Home Laboratory Computer Lab Theses

2 Optimization

”The great watershed in optimization is not between linearity and nonlinearity, but convexity and nonconvexity”

R. Tyrrell Rockafellar

The purpose of this chapter is to provide some background on the nonlinear programming (NLP) problem. Indeed, an exhaustive discussion of both theoretical and practical aspects of nonlinear programming can be the subject matter of an entire book.

There are several reasons for studying nonlinear programming in an optimal control class:

anyone interested in optimal control should know about a number of fundamental results in nonlinear programming; in fact, as optimal control problems are optimization problems in (infinite-dimensional) functional spaces while NLP problems are optimization problems in Euclidean spaces, optimal control can be seen as a generalization of nonlinear programming;
NLP techniques are used routinely and are particularly efficient in solving optimal control problems (in the case of a discrete control problem, i.e. when the controls are exerted at discrete points, the problem can be directly stated as a NLP problem; in a continuous control problem, on the other hand, i.e. when the controls are functions to be exerted over a prescribed planning horizon, an approximate solution can be found by solving a NLP problem).

2.1 Mathematical notation and useful definitions

We start by recalling some useful mathematical definition that will be used in the following.

Definition 2.1: Convex Set

A set $X \subset ℝ^{n}$ is said to be convex if $α x + (1 - α) y \in X \forall x, y \in X α \in [0, 1]$

Gemotrically convexity of a set means that every line segment joining two elements of the set belongs to the set

Definition 2.2: Convex Function

A function $f : X \subset ℝ^{n} \to ℝ$ is convex if $X$ is convex and holds that:

f (α x + (1 - α) y) \leq α f (x) + (1 - α) f (y) \forall x, y \in X α \in [0, 1]

Remark 2.1

A function is said strictly convex if the convexity condition holds with strict inequality, that is :

f (α x + (1 - α) y) < α f (x) + (1 - α) f (y) \forall x, y \in X α \in [0, 1]

Geometrically it means that the line segment joining $(x, f (x))$ and $(y, f (y))$ lies above f. Figures 2.1 , 2.2 illustrate the previous definitions.

Figure 2.2:: Illustration of convex function

Definition 2.3: Cone, Convex Cone

A nonempty set $C \subset ℝ^{n}$ is said to be a cone if, for every point $x \in C$ ,

α x \in C \forall α \geq 0

If, in addition, $C$ is convex then it is said to be a convex cone.

Definition 2.4: Positive definite matrix

A symmetric square matrix $A \in ℝ^{n}$ is said to be positive definite if $x ⊤ A x > 0 \forall x \in ℝ^{n}$

Remark 2.2

Equivalently, positive definiteness implies that all the eigenvalues are positive. In the same way we define semi-positive definite matrices. A symmetric square matrix $A \in ℝ^{n}$ is said to be positive definite if $x ⊤ A x \geq 0 \forall x \in ℝ^{n}$ . To denote positive and semi-positive definite matrices we use the symbol $A ≻ 0$ and $A ≽ 0$ respectively. If a symmetric matrix has both positive and negative eigenvalues is said to be indefinite.

Definition 2.5: small o

For a function $f : ℝ \to ℝ$ we write $f (x) = o (g (x))$ if and only if there exists a neighbourhood of $0$ (i.e. $ℬ_{𝜖} (0)$ ) and a function $c : ℬ_{𝜖} (0) \to ℝ$ with $\lim_{x \to 0} c (x) = 0$ such that:

| f (x) | \leq c (x) g (x) \forall x \in ℬ_{𝜖} (0)

Practically it means that if $f (x) = o (g (x))$ then $f$ shrinks faster than $g$ as $x$ goes to zero. For example $x^{2} = o (x)$ . Note that from the definition we also have that $\lim_{x \to 0} \frac{o (g (x))}{g (x)} = 0$

Definition 2.6

Matrix norm

2.2 Problem statement

We shall consider the following NLP problem:

Definition 2.7

Nonlinear Program

\begin{matrix} \min_{x} & f (x) \\ s.t. & g (x) \leq 0 \\ h (x) = 0 \\ x \in X \end{matrix}

where $X$ is a subset of $ℝ^{n_{x}}$ , $x$ is a vector of $n_{x}$ components $x_{1}, . . ., x_{n}$ , and $f : X \to ℝ$ , $g : X \to ℝ^{n_{g}}$ and $h : X \to ℝ^{n_{h}}$ are defined on $X$ .

Remark 2.3

The function $f$ is usually called the objective function or criterion function or performance index. Each of the constraints $g_{i} (x) \leq 0$ , $i = 1, . . ., n_{g}$ , is called an inequality constraint, and each of the constraints $h_{i} (x) = 0$ , $i = 1, . . ., n_{h}$ , is called an equality constraint.

Remark 2.4

Note also that the set $X$ typically includes lower and upper bounds on the variables. The reason for separating variable bounds from the other inequality constraints is that they can play a useful role in some algorithms, i.e. they are handled in a specific way.

Remark 2.5

A vector $x \in X$ satisfying all the constraints is called a feasible solution to the problem. The collection of all such points forms the feasible region. The NLP problem, then, is to find a feasible point $(x ⋆$ such that $f (x) \geq f ((x ⋆)$ for each feasible point $(x$ .

Needless to say, a NLP problem can be stated as a maximization problem, and the inequality constraints can be written in the form $g (x) \geq 0$ .

Example 2.1. Consider the following problem:

\begin{matrix} \min_{x} & {(x_{1} - 2)}^{2} + {(x_{2} - 1)}^{2} \\ s.t. & x_{1}^{2} - x_{2} \leq 0 \\ x_{1} + x_{2} - 2 \leq 0 \\ - x_{1} - 1 \leq 0 \end{matrix}

The objective function and the three inequality constraints are:

\begin{matrix} f (x_{1}, x_{2}) & = & {(x_{1} - 2)}^{2} + {(x_{2} - 1)}^{2} \\ g_{1} (x_{1}, x_{2}) & = & x_{1}^{2} - x_{2} \\ g_{2} (x_{1}, x_{2}) & = & x_{2} + x_{1} - 2 \\ g_{3} (x_{1}, x_{2}) & = & - x_{1} - 1 \end{matrix}

Figure 2.3 illustrates the feasible region. The problem, then, is to find a point in the feasible region with the smallest possible value of ${(x_{1} - 2)}^{2} + {(x_{2} - 1)}^{2}$ . Note that points $(x_{1}, x_{2})$ with ${(x_{1} - 2)}^{2} + {(x_{2} - 1)}^{2} = c$ are circles with radius $c$ and centre $(2, 1)$ . This circle is called the contour of the objective function having the value $c$ . In order to minimize $c$ , we must find the circle with the smallest radius that intersects the feasible region. As shown in figure 2.3 , the smallest circle corresponds to $c = 1$ and intersects the feasible region at the point $(1, 1)$ . Hence, the optimal solution occurs at the point $(1, 1)$ and has an objective value equal to $1$ .

The graphical approach used in the above example, i.e. find an optimal solution by determining the objective contour with the smallest objective value that intersects the feasible region, is only suitable for small problems. It becomes intractable for problems containing contours of the objective function in more than three variables as well as for problems having complicated objective and/or constraint functions.

pict — Figure 2.3:: Geometric solution of a nonlinear problem.

This chapter is organized as follows:

we start by defining what is meant by optimality and give conditions under which a minimum (or a maximum) exists for a nonlinear program;
then, we discuss the special properties of convex programs;
then, both necessary and sufficient conditions of optimality are presented for NLP problems;
successively, we consider unconstrained problems, problems with inequality constraints, and problems with both equality and inequality constraints;
finally, several numerical optimization techniques will be presented, which are instrumental to solve a great variety of NLP problems.

2.3 Definitions

A variety of different definitions of optimality are used in different contexts. It is important to understand fully each definition and the context within which it is appropriately used.

2.3.1 Infimum and Supremum

Let $S \subset ℝ$ be a nonempty set.

Definition 2.8: Infimum/Supremum

The infimum of $S$ , denoted as $\inf S$ , provided it exists, is the greatest lower bound for $S$ , i.e. a number $α$ satisfying:

$z \geq α \forall z \in S$ ,
$\forall \bar{α} > α, \exists z \in S$ such that $z < \bar{α}$ .

Similarly, the supremum of $S$ , denoted as $\sup S$ , provided it exists, is the least upper bound for $S$ , i.e. a number $α$ satisfying:

$z \leq α \forall z \in S$ ,
$\forall \bar{α} < α, \exists z \in S$ such that $z > \bar{α}$ .

The first question one may ask concerns the existence of infima and suprema in $ℝ$ . In particular, one cannot prove that in $ℝ$ every set bounded from above has a supremum and every set bounded from below has an infimum. This is an axiom, known as the axiom of completeness:

Axiom 2.1: of Completeness

If a nonempty subset of real numbers has an upper bound, then it has a least upper bound. If a nonempty subset of real numbers has a lower bound, it has a greatest lower bound.

Remark 2.6

It is important to note that the real number $\inf S$ (respectively $\sup S$ ), with $S$ a nonempty set in $ℝ$ bounded from below (respectively from above), although it exist, need not be an element of $S$ .

Example 2.2. Let $S = (0, + \infty) = {z \in ℝ : z > 0}$ . Clearly, $\inf S = 0$ and $0 \notin S$ .

Notation 2.1

Let $S : = f (x) : x \in D$ be the image of the feasible set $D \subset ℝ$ of an optimization problem under the objective function $f$ . Then, the notation

\begin{matrix} \inf_{x \in D} f (x) & or & \inf {f (x) : x \in D} \end{matrix}

refers to the number $\inf S$ . Likewise, the notation

\begin{matrix} \sup_{x \in D} f (x) & or & \sup {f (x) : x \in D} \end{matrix}

refers to $\sup S$ .

Clearly, the numbers $\inf S$ and $\sup S$ may not be attained by the value $f (x)$ at any $x \in D$ . This is illustrated in an example below.

Example 2.3. Clearly, $\inf {\exp (x) : x \in (0, + \infty)} = 1$ , but $\exp (x) > 1$ for all $x \in (0, + \infty)$ .

Remark 2.7

By convention, the infimum of an empty set is $+ \infty$ while the supremum of an empty set is $- \infty$ . That is, if the values $\pm \infty$ are allowed, then infima and suprema always exist.

2.3.2 Minimum and Maximum

Consider the standard problem formulation

\min_{x \in D} f (x)

where $D \subset ℝ^{n}$ denotes the feasible set. Any $x \in D$ is said to be a feasible point. Conversely, any $x \in ℝ^{n} ∖ D : = {x \in ℝ^{n} : x \notin D}$ is said to be an infeasible point.

Definition 2.9: (Global) Minimum, Strict (Global) Minimum

A point $x ⋆ \in D$ is said to be a (global) ^a minimum of $f$ on $D$ if

f (x) \geq f (x ⋆) \forall x \in D

(2.1)

i.e. a minimum is a feasible point whose objective function value is less than or equal to the objective function value of all other feasible points. It is said to be a strict (global) minimum of $f$ on $D$ if

f (x) > f (x ⋆) \forall x \in D, x \neq x ⋆

A (global) maximum is defined by reversing the inequality in the above definition:

Definition 2.10: (Global) Maximum, Strict (Global) Maximum

A point $x ⋆ \in D$ is said to be a (global) maximum of $f$ on $D$ if

f (x) \leq f (x ⋆) \forall x \in D

(2.2)

It is said to be a strict (global) maximum of $f$ on $D$ if

f (x) < f (x ⋆) \forall x \in D, x \neq x ⋆

Remark 2.8

The important distinction between minimum (respectively maximum) and infimum (respectively supremum) is that the value $\min {f (x) : x \in D}$ must be attained at one or more points $x \in D$ , whereas the value $\inf {f (x) : x \in D}$ does not necessarily have to be attained at any points $x \in D$ . Yet, if a minimum (respectively maximum) exists, then its optimal value will equal the infimum (respectively supremum).

Remark 2.9

Note also that if a minimum exists, it is not necessarily unique. That is, there may be multiple or even an infinite number of feasible points that satisfy the inequality 2.1 and are thus minima.

Notation 2.2

Since there is in general a set of points that are minima, the notation

\arg \min {f (x) : x \in D} : = {x \in D : f (x) = \inf {f (x) : x \in D}}

is introduced to denote the set of minima. This is a (possibly empty) set ^a ^a The notation $\bar{x} = \arg \min {f (x) : x \in D}$ is also used by some authors. In this case, $\arg \min {f (x) : x \in D}$ should be understood as a function returning a point $\bar{x}$ that minimizes $f$ on $D$ . in $ℝ^{n}$ .

Remark 2.10

A minimum $x ⋆$ is often referred to as an optimal solution, a global optimal solution, or simply a solution of the optimization problem. The real number $f (x ⋆)$ is known as the (global) optimal value or optimal solution value.

Notation 2.3

Regardless of the number of minima, there is always a unique real number that is the optimal value (if it exists). The notation $m i n {f (x) : x \in D}$ is used to refer to this real value.

Unless the objective function $f$ and the feasible set $D$ possess special properties (e.g. convexity), it is usually very hard to devise algorithms that are capable of locating or estimating a global minimum or a global maximum with certainty. This motivates the definition of local minima and maxima, which, by the nature of their definition in terms of local information, are much more convenient to locate with an algorithm.

Definition 2.11: Local Minimum, Strict Local Minimum

A point $x ⋆ \in D$ is said to be a local minimum of $f$ on $D$ if

\exists 𝜀 > 0 such that f (x) \geq f (x ⋆) \forall x \in ℬ_{𝜀} (x ⋆) \cap D

A point $x ⋆ \in D$ is said to be a strict local minimum of $f$ on $D$ if

\exists 𝜀 > 0 such that f (x) > f (x ⋆) \forall x \in ℬ_{𝜀} (x ⋆) \cap D

The qualifier ’local’ originates from the requirement that $x ⋆$ be a minimum only for those feasible points in a neighbourhood around the local minimum.

Remark 2.11

Trivially, the property of $x ⋆$ being a global minimum implies that $x ⋆$ is also a local minimum because a global minimum is local minimum with $𝜀$ set arbitrarily large.

Definition 2.12: Local Maximum, Strict Local Maximum

A point $x ⋆ \in D$ is said to be a local maximum of $f$ on $D$ if

\exists 𝜀 > 0 such that f (x) \leq f (x ⋆) \forall x \in ℬ_{𝜀} (x ⋆) \cap D

A point $x ⋆ \in D$ is said to be a strict local maximum of $f$ on $D$ if

\exists 𝜀 > 0 such that f (x) < f (x ⋆) \forall x \in ℬ_{𝜀} (x ⋆) \cap D

The qualifier ’local’ originates from the requirement that $x ⋆$ be a maximum only for those feasible points in a neighbourhood around the local maximum.

Remark 2.12

Trivially, the property of $x ⋆$ being a global maximum implies that $x ⋆$ is also a local maximum because a global maximum is local maximum with $𝜀$ set arbitrarily large.

Remark 2.13

It is important to note that the concept of a global minimum or a global maximum of a function on a set is defined without the notion of a distance (or a norm in the case of a vector space). In contrast, the definition of a local minimum or a local maximum requires that a distance be specified on the set of interest. In $ℝ_{x}^{n}$ , norms are equivalent, and it is readily shown that local minima (respectively maxima) in ( $ℝ_{x}^{n}$ , $∥ \cdot ∥_{α}$ ) match local minima (respectively maxima) in (( $ℝ_{x}^{n}$ , $∥ \cdot ∥_{β}$ ), for any two arbitrary norms $∥ \cdot ∥_{α}$ and $∥ \cdot ∥_{β}$ in $ℝ_{x}^{n}$ (e.g. the Euclidean norm $∥ \cdot ∥_{2}$ and the infinite norm $∥ \cdot ∥_{\infty}$ ). Yet, this nice property does not hold in linear functional spaces, as those encountered in problems of the calculus of variations and optimal control.

Figure ?? illustrates the various definitions of minima and maxima. Point $x_{1}$ is the unique global maximum; the objective value at this point is also the supremum. Points $a$ , $x_{2}$ , and $b$ are strict local minima because there exists a neighbourhood around each of these point for which $a$ , $x_{2}$ , or $b$ is the unique minimum (on the intersection of this neighbourhood with the feasible set $D$ ). Likewise, point $x_{3}$ is a strict local maximum. Point $x_{4}$ is the unique global minimum; the objective value at this point is also the infimum. Finally, point $x_{5}$ is simultaneously a local minimum and a local maximum because there are neighbourhoods for which the objective function remains constant over the entire neighbourhood; it is neither a strict local minimum, nor a strict local maximum.

Example 2.4. Consider the function

f (x) = {\begin{matrix} + 1 & if x < 0 \\ - 1 & otherwise \end{matrix}

The point $x^{⋆} = 0$ is a local minimum for

\min_{x \in [- 2, 2]} f (x)

with value $f (x^{⋆}) = - 1$ . The optimal value is $- 1$ and $\arg \min {f (x) : x \in [- 2, 2]} = [0, 2]$ .

2.3.3 Existence of Minima and Maxima

A crucial question when it comes to optimizing a function on a given set, is whether a minimizer or a maximizer exist for that function in that set. Strictly, a minimum or maximum should only be referred to when it is known to exist.

Figure ?? illustrates three instances where a minimum does not exist. In figure ??a, the infimum of $f$ over $S : = (a, b)$ is given by $f (b)$ , but since $S$ is not closed and, in particular, $b \notin S$ , a minimum does not exist. In figure ??b, the infimum of $f$ over $S : = [a, b]$ is given by the limit of $f (x)$ as $x$ approaches $c$ from the left, i.e. $\inf {f (x) : x \in S} = l i m_{x \to c^{-}} f (x)$ . However, because $f$ is discontinuous at $c$ , a minimizing solution does not exist. Finally, figure ??c illustrates a situation within which $f$ is unbounded over the unbounded set $S : = {x \in ℝ : x \geq a}$ .

We now formally state and prove the result that if $S$ is nonempty, closed, and bounded, and if $f$ is continuous on $S$ , then, unlike the various situations of figure ??, a minimum exists.

Theorem 2.1: Weierstrass’ Theorem

Let $S$ be a nonempty compact set and let $f : S \to ℝ$ be continuous on $S$ . Then, the problem $\min {f (x) : x \in S}$ attains its minimum, that is, there exists a minimizing solution to this problem.

Proof. Since $f$ is continuous on $S$ and $S$ is both closed and bounded, $f$ is bounded below on $S$ . Consequently, since $S \neq \emptyset$ , there exists a greatest lower bound $α : = \inf {f (x) : x \in S}$ (see previous axiom). Now, let $0 < 𝜀 < 1$ , and consider the set $S_{k} : = {x \in S : α \leq f (x) \leq α + 𝜀^{k}}$ for $k = 1, 2, . . .$ . By the definition of an infimum, $S_{k} \neq \emptyset$ for each $k$ , and so we may construct a sequence of points ${x_{k}} \subset S$ by selecting a point $x_{k}$ for each $k = 1, 2, . . .$ . Since $S$ is bounded, there exists a convergent subsequence ${x_{k}}_{𝒦} \subset S$ indexed by the set $𝒦 \subset ℕ$ . Let $\bar{x}$ denote its limit. By the closedness of $S$ , we have $\bar{x} \in S$ and by the continuity of $f$ on $S$ , since $α \leq f (x k) \leq α + 𝜀^{k}$ , we have $α = l i m_{k \to \infty, k \in 𝒦} f (x k) = f (\bar{x})$ . Hence, we have shown that there exist a solution $\bar{x} \in S$ such that $f (\bar{x}) = α = \inf {f (x) : x \in S}$ , i.e. $\bar{x}$ is a minimizing solution. □

Remark 2.14

The hypotheses of theorem 2.1 can be justified as follows:

the feasible set must be nonempty, otherwise there are no feasible points at which to attain the minimum;
the feasible set must contain its boundary points, which is ensured by assuming that the feasible set is closed;
the objective function must be continuous on the feasible set, otherwise the limit at a point may not exist or be different from the value of the function at that point;
the feasible set must be bounded because otherwise even a continuous function can be unbounded on the feasible set.

Example 2.5. Theorem 2.1 establishes that a minimum (and a maximum) of

\min_{x \in [- 1, 1]} x^{2}

exists, since $[- 1, 1]$ is a nonempty compact set and $x \mapsto x^{2}$ is a continuous function on $[- 1, 1]$ . On the other hand, minima can still exist even though the set is not compact or the function is not continuous, for theorem 2.1 only provides a sufficient condition. This is the case for the problem

\min_{x \in (- 1, 1)} x^{2}

which has a minimum at $x = 0$ .

Example 2.6. Consider the NLP problem

\begin{matrix} \min_{x} & {(x_{1} - 3)}^{2} + {(x_{2} - 2)}^{2} \\ s.t. & x_{1}^{2} - x_{2} - 3 \leq 0 \\ x_{2} - 1 \leq 0 \\ - x_{1} \leq 0 \end{matrix}

The objective function being continuous and the feasible region being nonempty, closed and bounded, the existence of a minimum to this problem directly follows from theorem 2.1 .

2.4 Convex programming

A particular class of nonlinear programs is that of convex programs:

Definition 2.13: Convex Program

Let $C$ be a nonempty convex set in $ℝ^{n}$ , and let $f : C \to ℝ$ be convex on $C$ . Then,

\min_{x \in C} f (x)

is said to be a convex program (or a convex optimization problem).

Convex programs possess nicer theoretical properties than general nonconvex problems. Theorem 2.2 is a fundamental result in convex programming:

Theorem 2.2

Let $x ⋆$ be a local minimum of a convex program. Then, $x ⋆$ is also a global minimum.

Proof. $x ⋆$ being a local minimum,

\exists 𝜀 > 0 such that f (x) \geq f (x ⋆), \forall x \in ℬ_{𝜀} (x ⋆)

By contradiction, suppose that $x ⋆$ is not a global minimum. Then,

\exists \bar{x} \in C such that f (\bar{x}) < f (x ⋆)

(2.3)

Let $λ \in (0, 1)$ be chosen such that $y : = λ \bar{x} + (1 - λ) x ⋆ \in ℬ_{𝜀} (x ⋆)$ . By convexity of $C$ , $y$ is in $C$ . Next, by convexity of $f$ on $C$ and 2.3,

f (y) \leq λ f (\bar{x}) + (1 - λ) f (x ⋆) < λ f (x ⋆) + (1 - λ) f (x ⋆) = f (x ⋆)

hence contradicting the assumption that $x ⋆$ is a local minimum. □

Example 2.7. Consider once again the NLP problem

\begin{matrix} \min_{x} {(x_{1} - 3)}^{2} + {(x_{2} - 2)}^{2} \\ s.t. x_{1}^{2} - x_{2} - 3 \leq & 0 \\ x_{2} - 1 \leq & 0 \\ - x_{1} & \leq 0 \end{matrix}

The objective function $f$ and the inequality constraints $g_{1}$ , $g_{2}$ and $g_{3}$ being convex, every local solution to this problem is also a global solution by theorem 2.2 . Henceforth, $x ⋆ = (1, 2)$ is a global solution and the global solution value is $f (x ⋆) = 4$ .

In convex programming, any local minimum is therefore a global minimum. This is a powerful result that makes any local optimization algorithm a global optimization algorithm when applied to a convex optimization problem. Yet, theorem 2.2 only gives a sufficient condition for that property to hold. That is, a nonlinear program with nonconvex participating functions may not necessarily have local minima that are not global minima.

2.5 Unconstrained problems

Definition 2.14: Unconstrained Programs

An unconstrained program is a problem of the form to minimize (or maximize) $f (x)$ without any constraints on the variables $x$ :

\min {f (x) : x \in ℝ^{n_{x}}}

Remark 2.15

Note that, being the feasible domain of $x$ unbounded, theorem 2.1 does not apply. Thus, one does not know with certainty whether a minimum actually exists for that problem ^a ^a For unconstrained optimization problems, the existence of a minimum can actually be guaranteed if the objective function is such that $\lim_{∥ x ∥ \to + \infty} f (x) = + \infty$ (O-coercive function). . Moreover, even if the objective function is convex, one such minimum may not exist (think of $f : x \mapsto \exp x$ !). Hence, we shall proceed with the theoretically unattractive task of seeking minima and maxima of functions which need not have them!

Given a point $x \in ℝ^{n_{x}}$ , necessary conditions help determine whether or not a point is a local or a global minimum of a function $f$ . For this purpose, we are mostly interested in obtaining conditions that can be checked algebraically.

Definition 2.15: Descent Direction

Suppose that $f : ℝ^{n_{x}} \to ℝ$ is continuous at $\bar{x}$ . A vector $d \in ℝ^{n_{x}}$ is said to be a descent direction, or an improving direction, for $f$ at $\bar{x}$ if

\exists δ > 0 : f (\bar{x} + λ d) < f (\bar{x}) \forall λ \in (0, δ)

Moreover, the cone of descent directions at $\bar{x}$ , denoted by $ℱ (\bar{x})$ , is given by

ℱ (\bar{x}) : = {d : \exists δ > 0 such that f (\bar{x} + λ d) < f (\bar{x}) \forall λ \in (0, δ)}

The definition 2.15 provides a geometrical characterization for a descent direction. Yet, an algebraic characterization for a descent direction would be more useful from a practical point of view. In response to this, let us assume that $f$ is differentiable and define the following set at $\bar{x}$ :

ℱ_{0} (\bar{x}) : = {d : \nabla f {(\bar{x})}^{⊤} d < 0}

This is illustrated in figure ?? where the half-space $ℱ_{0} (\bar{x})$ and the gradient $\nabla f (\bar{x})$ are translated from the origin to $\bar{x}$ for convenience.

Lemma 2.1 proves that every element $d \in ℱ_{0} (\bar{x})$ is a descent direction at $\bar{x}$ .

Lemma 2.1: Algebraic Characterization of a Descent Direction

Suppose that $f : ℝ^{n_{x}} \to ℝ$ is differentiable at $\bar{x}$ . If there exists a vector $d$ such that $\nabla f {(\bar{x})}^{⊤} d < 0$ , then $d$ is a descent direction for $f$ at $\bar{x}$ . That is,

ℱ_{0} (\bar{x}) \subseteq ℱ (\bar{x})

Proof. $f$ being differentiable at $\bar{x}$ ,

f (\bar{x} + λ d) = f (\bar{x}) + λ \nabla f {(\bar{x})}^{⊤} d + o (λ)

where $\lim_{λ \to 0} \frac{o (λ)}{λ} = 0$ . Rearranging the terms and dividing by $λ \neq 0$ , we get

\frac{f (\bar{x} + λ d) - f (\bar{x})}{λ} = \nabla f {(\bar{x})}^{⊤} d + \frac{o (λ)}{λ}

Since $\nabla f {(\bar{x})}^{⊤} d < 0$ and $\lim_{λ \to 0} \frac{o (λ)}{λ} = 0$ , there exists a $δ > 0$ such that $\nabla f {(\bar{x})}^{⊤} d + \frac{o (λ)}{λ} < 0$ for all $λ \in (0, δ)$ . □

We are now ready to derive a number of necessary conditions for a point to be a local minimum of an unconstrained optimization problem.

Theorem 2.3: First-Order Necessary Condition for a Local Minimum

Suppose that $f : ℝ^{n_{x}} \to ℝ$ is differentiable at $x ⋆$ . If $x ⋆$ is a local minimum, then $\nabla f (x ⋆) = 0$ .

Proof. The proof proceeds by contraposition. Suppose that $\nabla f (x ⋆) \neq 0$ . Then, letting $d = - \nabla f (x ⋆)$ , we get $\nabla f (x ⋆) ⊤ d = - ∥ \nabla f (x ⋆) ∥ 2 < 0$ . By lemma 2.1,

\exists δ > 0 : f (x ⋆ + λ d) < f (x ⋆) \forall λ \in (0, δ)

hence contradicting the assumption that $x ⋆$ is a local minimum for $f$ . □

Remark 2.16: Obtaining Candidate Solution Points

The above condition is called a first-order necessary condition because it uses the first-order derivatives of $f$ . This condition indicates that the candidate solutions to an unconstrained optimization problem can be found by solving a system of $n_{x}$ algebraic (nonlinear) equations. Points $\bar{x}$ such that $\nabla f (\bar{x}) = 0$ are known as stationary points. Yet, a stationary point need not be a local minimum; it could very well be a local maximum or even a saddle point.

Example 2.8. Consider the problem

\min_{x \in ℝ} - x^{6} + x^{4} + x^{3} - x^{2}

The gradient vector of the objective function is given by

\nabla f (x) = - 6 x^{5} + 4 x^{3} + 3 x^{2} - 2 x

which has three distinct roots $x_{1}^{⋆}$ , $x_{2}^{⋆}$ and $x_{3}^{⋆}$ . Out of these values, $x_{2}^{⋆}$ gives the smallest cost value. Yet, we cannot declare $x_{2}^{⋆}$ to be the global minimum because we do not know whether a (global) minimum exists for this problem. Indeed, as shown in figure 2.4 , none of the stationary points is a global minimum because $f$ decreases to $- \infty$ as $| x | \to \infty$ .

More restrictive necessary conditions can also be derived in terms of the Hessian matrix $\nabla^{2} f$ whose elements are the second-order derivatives of $f$ . One such second-order condition is given below.

Theorem 2.4: Second-Order Necessary Conditions for a Local Minimum

Suppose that $f : ℝ^{n_{x}} \to ℝ$ is twice differentiable at $x ⋆$ . If $x ⋆$ is a local minimum, then $\nabla f (x ⋆) = 0$ and $\nabla^{2} f (x ⋆)$ is positive semidefinite.

Proof. Consider an arbitrary direction $d$ . Then, from the differentiability of $f$ at $x ⋆$ , we have

f (x ⋆ + λ d) = f (x ⋆) + λ \nabla f (x ⋆) ⊤ d + \frac{λ^{2}}{2} d ⊤ \nabla^{2} f (x ⋆) d + o (λ^{2})

where $\lim_{λ \to 0} \frac{o (λ^{2})}{λ^{2}} = 0$ . Since $x ⋆$ is a local minimum, from Theorem 2.3, $\nabla f (x ⋆) = 0$ . Rearranging the terms and dividing by $λ^{2}$ , we get

\frac{f (x ⋆ + λ d)}{-} = \frac{1}{2} d ⊤ \nabla^{2} f (x ⋆) d + \frac{o (λ^{2})}{λ^{2}}

Since $x ⋆$ is a local minimum, $f (x ⋆ + λ d) \geq f (x ⋆)$ for $λ$ sufficiently small. By taking the limit as $λ \to 0$ , it follows that $d ⊤ \nabla^{2} f (x ⋆) d \geq 0$ . Since $d$ is arbitrary, $\nabla^{2} f (x ⋆)$ is therefore positive semidefinite. □

Example 2.9. Consider the problem

\min_{x \in ℝ^{2}} \frac{1}{2} x ⊤ [\begin{matrix} 3 & - 2 \\ - 2 & - 1.5 \end{matrix}] x = \frac{1}{2} x ⊤ A x

The gradient vector of the objective function is given by

\nabla f (x) = A x

so that the only stationary point in $ℝ^{2}$ is $\bar{x} = (0, 0)$ . Now, consider the Hessian matrix of the objective function at $\bar{x}$ :

\nabla^{2} f (\bar{x}) = A = (\begin{matrix} 3 & - 2 \\ - 2 & - 1.5 \end{matrix}) \forall x \in ℝ^{2}

It is easily checked that $\nabla^{2} f (\bar{x})$ is indefinite, $eig (A) = {3.8, - 2.3}$ . Therefore, by Theorem 2.4 , the stationary point $\bar{x}$ is not a (local) minimum (nor is it a local maximum). Such stationary points are called saddle points (see Figure 2.5 ).

The conditions presented in theorems 2.3 and 2.4 are necessary conditions. That is, they must hold true at every local optimal solution. Yet, a point satisfying these conditions need not be a local minimum. Theorem 2.5 gives sufficient conditions for a stationary point to be a global minimum point, provided the objective function is convex on $ℝ^{n_{x}}$ .

Theorem 2.5: First-Order Sufficient Conditions for a Global Minimum for Convex Functions

Suppose that $f : ℝ^{n_{x}} \to ℝ$ is differentiable at $x ⋆$ and convex on $ℝ^{n_{x}}$ . If $\nabla f (x ⋆) = 0$ , then $x ⋆$ is a global minimum of $f$ on $ℝ^{n_{x}}$ .

Proof. $f$ being convex on $ℝ^{n_{x}}$ and differentiable at $x ⋆$ , by Theorem 2.6, we have

f (x) \geq f (x ⋆) + \nabla f (x ⋆) ⊤ [x - x ⋆] \forall x \in ℝ^{n_{x}}

But since $x ⋆$ is a stationary point,

f (x) \geq f (x ⋆) \forall x \in ℝ^{n_{x}}

□

Theorem 2.6: First-Order Condition of Convexity

Let $C$ be a convex set in $ℝ^{n_{x}}$ with a nonempty interior, and let $f : C \to ℝ$ be a function. Suppose $f$ is continuous on $C$ and differentiable on $int (C)$ . Then $f$ is convex on $int (C)$ if and only if

f (y) \geq f (x) + \nabla f {(x)}^{⊤} [y - x]

holds for any two points $x, y \in C$ .

Proof. ... □

Theorem 2.7: Second-Order Condition of Convexity

Let $C$ be a convex set in $ℝ^{n_{x}}$ with a nonempty interior, and let $f : C \to ℝ$ be a function. Suppose $f$ is continuous on $C$ and twice differentiable on $int (C)$ . Then $f$ is convex on $int (C)$ if and only if

\nabla^{2} f (x) ≽ 0 \forall x \in C

Proof. ... □

Remark 2.17

The first-order condition of strict convexity is:

f (y) > f (x) + \nabla f {(x)}^{⊤} [y - x] \forall x, y \in C

While the second-order condition is:

\nabla^{2} f (x) ≻ 0 \forall x \in C

The convexity condition required by theorem 2.5 is actually very restrictive, in the sense that many practical problems are nonconvex. In theorem 2.8, we give sufficient conditions for characterizing a local minimum point, provided the objective function is strictly convex in a neighborhood of that point.

Theorem 2.8: Second-Order Sufficient Conditions for a Strict Local Minimum

Suppose that $f : ℝ^{n_{x}} \to ℝ$ is twice differentiable at $x ⋆$ . If $\nabla f (x ⋆) = 0$ and $\nabla^{2} f (x ⋆)$ is positive definite, then $x ⋆$ is a local minimum of $f$ .

Proof. $f$ being twice differentiable at $x ⋆$ , we have

f (x ⋆ + d) = f (x ⋆) + \nabla f (x ⋆) ⊤ d + \frac{1}{2} d ⊤ \nabla^{2} f (x ⋆) d + o (∥ d ∥^{2})

for each $d \in ℝ^{n_{x}}$ , where $\lim_{∥ d ∥ \to 0} \frac{o (∥ d ∥^{2})}{∥ d ∥^{2}} = 0$ . Let $λ^{L}$ denote the smallest eigenvalue of $\nabla^{2} f (x ⋆)$ . Then, $\nabla^{2} f (x ⋆)$ being positive definite we have $λ^{L} > 0$ and $d ⊤ \nabla^{2} f (x ⋆) d \geq λ^{L} ∥ d ∥^{2}$ . Moreover, from $\nabla f (x ⋆) = 0$ we get

f (x ⋆ + d) - f (x ⋆) \geq [\frac{λ^{L}}{2} + \frac{o (∥ d ∥^{2})}{∥ d ∥^{2}}] ∥ d ∥^{2}

Since $\lim_{∥ d ∥ \to 0} \frac{o (∥ d ∥^{2})}{∥ d ∥^{2}} = 0$ ,

\exists η > 0 such that | \frac{o (∥ d ∥^{2})}{∥ d ∥^{2}} | < \frac{λ^{L}}{2} \forall d \in 𝔹_{η} (0)

and finally,

f (x ⋆ + d) - f (x ⋆) \geq \frac{λ^{L}}{2} ∥ d ∥^{2} > 0 \forall d \in 𝔹_{η} (0) ∖ {0}

i.e., $x ⋆$ is a strict local minimum of $f$ . □

Example 2.10. Consider the problem

\min_{x \in ℝ^{2}} \frac{1}{2} {[\begin{matrix} x_{1} - 1 \\ x_{2} - 2 \end{matrix}]}^{⊤} [\begin{matrix} 5 & - 2 \\ - 2 & 2 \end{matrix}] [\begin{matrix} x_{1} - 1 \\ x_{2} - 2 \end{matrix}] = \frac{1}{2} {[\begin{matrix} x_{1} - 1 \\ x_{2} - 2 \end{matrix}]}^{⊤} A [\begin{matrix} x_{1} - 1 \\ x_{2} - 2 \end{matrix}]

The gradient vector and Hessian matrix at $\bar{x} = (1, 2)$ are given by

\nabla f (\bar{x}) = A [\begin{matrix} x_{1} - 1 \\ x_{2} - 2 \end{matrix}] = 0

\nabla^{2} f (\bar{x}) = A ≻ 0

Note that $A ≻ 0$ since $eig (A) = {6, 1}$ and also that $\nabla^{2} f ≽ 0$ is the second-order characterization of strict convexity. Hence, by Theorem 2.8 , $\bar{x}$ is a local minimum of $f$ ( $\bar{x}$ is also a global minimum of $f$ on $ℝ^{2}$ since $f$ is convex). The objective function is pictured in Figure 2.6 below.

We close this subsection by re-emphasizing the fact that every local minimum of an unconstrained optimization problem $\min {f (x) : x \in ℝ^{n_{x}}}$ is a global minimum if $f$ is a convex function on $ℝ^{n_{x}}$ . Yet, convexity of $f$ is not a necessary condition for each local minimum to be a global minimum.

2.6 Problems with inequality constraints

In practice, few problems can be formulated as unconstrained programs. This is because the feasible region is generally restricted by imposing constraints on the optimization variables. In this section, we first present theoretical results for the problem that is to minimize $f (x)$ , subject to $x \in S$ for a general set $S$ (geometric optimality conditions). Then, we let $S$ be more specifically defined as :

Definition 2.16: Inequality Constrained Program

\begin{matrix} minimize & f (x) \\ subject to & g (x) \leq 0 \end{matrix}

The feasible region of the NLP is defined by a set of nonlinear inequalities $g_{i} (x) \leq 0$ , for the Inequality Constrained Program we derive the Karush-Kuhn-Tucker (KKT) conditions of optimality in the following.

2.6.1 Geometric Optimality Conditions

Definition 2.17: Feasible Direction

Let $S$ be a nonempty set in $ℝ^{n_{x}}$ . A vector $d \in ℝ^{n_{x}}$ , $d \neq 0$ , is said to be a feasible direction at $\bar{x} \in c l (S)$ ^a ^a cl(S) indicates the closure of the set S, that is the union of S with all of its limit points if

\exists δ > 0 such that \bar{x} + η d \in S \forall η \in (0, δ)

Moreover, the cone of feasible directions at $\bar{x}$ , denoted by $𝒟 (\bar{x})$ , is given by

𝒟 (\bar{x}) : = {d \neq 0 : \exists δ > 0 such that \bar{x} + η d \in S \forall η \in (0, δ)}

From the Definition 2.17 and from Lemma 2.1, it is clear that a small movement from $\bar{x}$ along a direction $d \in 𝒟 (\bar{x})$ leads to feasible points, whereas a similar movement along a direction $d \in ℱ_{0} (\bar{x})$ leads to solutions of improving objective value (see the definition 2.15). As shown in theorem ??, a (geometric) necessary condition for local optimality is that: ”Every improving direction is not a feasible direction”. This fact is illustrated in figure ?? where both the half-space $ℱ_{0} (\bar{x})$ and the cone $𝒟 (\bar{x})$ (see Definition 2.3) are translated from the origin to $\bar{x}$ for clarity.

Theorem 2.9: Geometric Necessary Condition for a Local Minimum

Let $S$ be a nonempty set in $ℝ^{n_{x}}$ and let $f : ℝ^{n_{x}} \to ℝ$ be a differentiable function. Suppose that $\bar{x}$ is a local minimizer of the problem to minimize $f (x)$ subject to $x \in S$ . Then, $ℱ_{0} (\bar{x}) \cap 𝒟 (\bar{x}) = \emptyset$ .

Proof. By contradiction, suppose that there exists a vector $d \in ℱ_{0} (\bar{x}) \cap 𝒟 (\bar{x})$ , $d \neq 0$ . Then, by lemma 2.1,

\exists δ_{1} > 0 such that f (\bar{x} + η d) < f (\bar{x}) \forall η \in (0, δ_{1})

Moreover, by the definition 2.17,

\exists δ_{2} > 0 such that \bar{x} + η d \in S \forall η \in (0, δ_{2})

Hence,

\exists x \in ℬ_{η} (\bar{x}) \cap S such that f (\bar{x} + η d) < f (\bar{x})

for every $η \in (0, \min {δ_{1}, δ_{2}})$ , which contradicts the assumption that $\bar{x}$ is a local minimum of $f$ on $S$ (see the definition 2.11). □

2.6.2 KKT Conditions

We now specify the feasible region as

S : = {x : g_{i} (x) \leq 0 \forall i = 1, \dots, n_{g}}

where $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ , are continuous functions. In the geometric optimality condition given by Theorem 2.12, $𝒟 (\bar{x})$ is the cone of feasible directions. From a practical viewpoint, it is desirable to convert this geometric condition into a more usable condition involving algebraic equations. As Lemma 2.2 indicates, we can define a cone $𝒟_{0} (\bar{x})$ in terms of the gradients of the active constraints at $\bar{x}$ , such that $𝒟_{0} (\bar{x}) \subseteq 𝒟 (\bar{x})$ . For this, we need the following:

Definition 2.18: Active Constraint, Active Set

Let $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ , and consider the set
$S : = {x : g_{i} (x) \leq 0, i = 1, \dots, n_{g}}$ . Let $\bar{x} \in S$ be a feasible point. For each $i = 1, \dots, n_{g}$ , the constraint $g_{i}$ is said to be active or binding at $\bar{x}$ if $g_{i} (\bar{x}) = 0$ ; it is said to be inactive at $\bar{x}$ if $g_{i} (\bar{x}) < 0$ . Moreover,

𝒜 (\bar{x}) : = {i : g_{i} (\bar{x}) = 0}

denotes the set of active constraints at $\bar{x}$ .

Lemma 2.2: Algebraic Characterization of a Feasible Direction

Let $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ be differentiable functions, and consider the set $S : = {x : g_{i} (x) \leq 0, i = 1, \dots, n_{g}}$ . For any feasible point $\bar{x} \in S$ , we have

𝒟_{0} (\bar{x}) : = {d : \nabla g_{i} {(\bar{x})}^{⊤} d < 0 \forall i \in 𝒜 (\bar{x})} \subseteq 𝒟 (\bar{x})

Proof. Suppose $𝒟_{0} (\bar{x})$ is nonempty, and let $d \in 𝒟_{0} (\bar{x})$ . Since $\nabla g_{i}^{⊤} (\bar{x}) d < 0$ for each $i \in 𝒜 (\bar{x})$ , then by Lemma 2.3, $d$ is a descent direction for $g_{i}$ at $\bar{x}$ , i.e.,

\exists δ_{2} > 0 such that g_{i} (\bar{x} + η d) < g_{i} (\bar{x}) = 0 \forall η \in (0, δ_{2}), \forall i \in 𝒜 (\bar{x})

Furthermore, since $g_{i} (\bar{x}) < 0$ and $g_{i}$ is continuous at $\bar{x}$ (since it is differentiable) for each $i \notin 𝒜 (\bar{x})$ ,

\exists δ_{1} > 0 such that g_{i} (\bar{x} + η d) < 0 \forall η \in (0, δ_{1}), \forall i \notin 𝒜 (\bar{x})

Furthermore, overall, it is clear that the points $\bar{x} + η d$ are in $S$ for all $η \in (0, \min {δ_{1}, δ_{2}})$ . Hence, by Definition 2.17, $d \in 𝒟 (\bar{x})$ . □

Remark 2.18

This Lemma together with Theorem 2.12 directly leads to the result that $ℱ_{0} (\bar{x}) \cap 𝒟_{0} (\bar{x}) = \emptyset$ for any local solution point $\bar{x}$ , i.e.,

\arg \min {f (x) : x \in S} \subset {x \in ℝ^{n_{x}} : ℱ_{0} (x) \cap 𝒟_{0} (x) = \emptyset}

The foregoing geometric characterization of local solution points applies equally well to either interior points
$int (S) : = {x \in ℝ^{n_{x}} : g_{i} (x) < 0, \forall i = 1, \dots, n_{g}}$ or boundary points being at the boundary of the feasible domain. At an interior point, in particular, any direction is feasible, and the necessary condition $ℱ_{0} (\bar{x}) \cap 𝒟_{0} (\bar{x}) = \emptyset$ reduces to $\nabla f (\bar{x}) = 0$ , which gives the same condition as in unconstrained optimization (see Theorem 2.3).

Note also that there are several cases where the condition $ℱ_{0} (\bar{x}) \cap 𝒟_{0} (\bar{x}) = \emptyset$ is satisfied by non-optimal points. In other words, this condition is necessary but not sufficient for a point $\bar{x}$ to be a local minimum of $f$ on $S$ . For instance, any point $\bar{x}$ with $\nabla g_{i} (\bar{x}) = 0$ for some $i \in 𝒜 (\bar{x})$ trivially satisfies the condition $ℱ_{0} (\bar{x}) \cap 𝒟_{0} (\bar{x}) = \emptyset$ . Another example is given below.

Example 2.11. Consider the problem

\begin{aligned} \min_{x \in ℝ^{2}} f (x) & : = x_{1}^{2} + x_{2}^{2} \\ s.t. & g_{1} (x) : = x_{1} \leq 0 \\ g_{2} (x) : = - x_{1} \leq 0 \end{aligned}

Clearly, this problem is convex and $x ⋆ = {(0, 0)}^{⊤}$ is the unique global minimum.

Now, let $\bar{x}$ be any point on the line $𝒞 : = {x : x_{1} = 0}$ . Both constraints $g_{1}$ and $g_{2}$ are active at $\bar{x}$ , and we have $\nabla g_{1} (\bar{x}) = - \nabla g_{2} (\bar{x}) = {(1, 0)}^{⊤}$ . Therefore, no direction $d \neq 0$ can be found such that $\nabla g_{1} {(\bar{x})}^{⊤} d < 0$ and $\nabla g_{2} {(\bar{x})}^{⊤} d < 0$ simultaneously, i.e., $𝒟_{0} (\bar{x}) = \emptyset$ . In turn, this implies that $ℱ_{0} (\bar{x}) \cap 𝒟_{0} (\bar{x}) = \emptyset$ is trivially satisfied for any point on $𝒞$ .

On the other hand, observe that the condition $ℱ_{0} (\bar{x}) \cap 𝒟 (\bar{x}) = \emptyset$ in Theorem 2.9 excludes all the points on $𝒞$ , but the origin, since a feasible direction at $\bar{x}$ is given, e.g., by $d = {(0, 1)}^{⊤}$ .

We now reduce the geometric necessary optimality condition $ℱ_{0} (\bar{x}) \cap 𝒟_{0} (\bar{x}) = \emptyset$ to a statement in terms of the gradients of the objective function and of the active constraints. The resulting first-order optimality conditions are known as the Karush-Kuhn-Tucker (KKT) necessary conditions. Beforehand, we introduce the important concepts of a regular point and of a KKT point.

Definition 2.19: Regular Point (for a Set of Inequality Constraints)

Let $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ be differentiable functions on $ℝ^{n_{x}}$ and consider the set $S : = {x \in ℝ^{n_{x}} : g_{i} (x) \leq 0, i = 1, \dots, n_{g}}$ . A point $\bar{x} \in S$ is said to be a regular point if the gradient vectors $\nabla g_{i} (\bar{x}), i \in 𝒜 (\bar{x})$ are linearly independent:

rank (\nabla g_{i} (\bar{x}), i \in 𝒜 (\bar{x})) = | 𝒜 (\bar{x}) |

Definition 2.20: KKT Point

Let $f : ℝ^{n_{x}} \to ℝ$ and $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ be differentiable functions. Consider the problem to minimize $f (x)$ subject to $g_{i} (x) \leq 0, i = 1, \dots, n_{g}$ . If a point $(\bar{x}, \bar{ν}) \in ℝ^{n_{x}} \times ℝ^{n_{g}}$ satisfies the conditions:

\begin{align} \nabla f (\bar{x}) + \nabla g (\bar{x}) \bar{ν} & = \nabla f (\bar{x}) + \sum_{i = 1}^{n_{g}} \nabla g_{i} (\bar{x}) {\bar{ν}}_{i} = 0 & (2.4a) \\ \bar{ν} & \geq 0 & (2.4b) \\ g (\bar{x}) & \leq 0 & (2.4c) \\ {\bar{ν}}^{⊤} g (\bar{x}) & = 0 & (2.4d) \end{align}

then $(\bar{x}, \bar{ν})$ is said to be a KKT point ^a ^a Note that $\nabla g = {\frac{\partial g}{\partial x}}^{⊤} = [\nabla g_{1} \dots \nabla g_{n_{g}}]$ .

Remark 2.19

The scalars $ν_{i}, i = 1, \dots, n_{g}$ are called the Lagrange multipliers. The condition 2.4a is called stationarity condition. The condition 2.4b is referred as dual feasibility (DF) condition; the condition 2.4c, i.e., the requirement that $\bar{x}$ be feasible, is called the primal feasibility (PF) condition; finally, the condition 2.4d is called the complementarity slackness (CS) condition ^a

{\bar{ν}}_{i} g_{i} (\bar{x}) = 0 for i = 1, \dots, n_{g}

Theorem 2.10: KKT Necessary Conditions

Let $f : ℝ^{n_{x}} \to ℝ$ and $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ be differentiable functions. Consider the problem to minimize $f (x)$ subject to $g_{i} (x) \leq 0, i = 1, \dots, n_{g}$ . If $x ⋆$ is a local minimum and a regular point of the constraints, then there exists a unique vector $ν ⋆$ such that $(x ⋆, ν ⋆)$ is a KKT point.

Proof. Since $x ⋆$ solves the problem, then there exists no direction $d \in ℝ^{n_{x}}$ such that $\nabla f {(x}^{⋆}) ⊤ d < 0$ and $\nabla g_{i} {(x}^{⋆}) ⊤ d < 0$ , $\forall i \in 𝒜 (x ⋆)$ simultaneously (see Remark 2.18). Let $A \in ℝ^{(| 𝒜 (x ⋆}) | + 1) \times n_{x}$ be the matrix whose rows are $\nabla f {(\bar{x})}^{⊤}$ and $\nabla g_{i} {(\bar{x})}^{⊤}$ , $i \in 𝒜 (x ⋆)$ . Clearly, the statement ${\exists d \in ℝ^{n_{x}} : A d < 0}$ is false and, by Corollary 2.1, there exists a nonzero vector $p \geq 0$ in $ℝ^{| 𝒜 (x ⋆}) | + 1$ such that $A^{⊤} p = 0$ . Denoting the components of $p$ by $u_{0}$ and $u_{i}$ for $i \in 𝒜 (x ⋆)$ , we get:

u_{0} \nabla f (x ⋆) + \sum_{i \in 𝒜 (x ⋆)} u_{i} \nabla g_{i} (x ⋆) = 0

where $u_{0} \geq 0$ and $u_{i} \geq 0$ for $i \in 𝒜 (x ⋆)$ and $(u_{0}, u 𝒜 (x ⋆)) \neq (0, 0)$ (here $u 𝒜 (x ⋆)$ is the vector whose components are the $u_{i}$ ’s for $i \in 𝒜 (x ⋆)$ ). Letting $u_{i} = 0$ for $i \notin 𝒜 (x ⋆)$ , we then get the conditions:

\begin{aligned} u_{0} \nabla f (x ⋆) + \nabla g (x ⋆) u & = u_{0} \nabla f (x ⋆) + \sum_{i = 1}^{n_{g}} u_{i} \nabla g_{i} (x ⋆) = 0 \\ u ⊤ g (x ⋆) & = 0 \\ (u_{0}, u) & \geq (0, 0) \\ (u_{0}, u) & \neq (0, 0) \end{aligned}

where $u$ is the vector whose components are $u_{i}$ for $i = 1, \dots, n_{g}$ . Note that $u_{0} \neq 0$ , since otherwise the assumption of linear independence of the active constraints at $x ⋆$ would be violated ^a ^a Note that if $u_{0} = 0$ then $u \neq 0$ and by Corollary 2.1 we have:

\sum_{i \in 𝒜 (x ⋆)} u_{i} \nabla g_{i} (x ⋆) = 0

with at least one $u_{i} \neq 0$ thus violating the linear independence assumption. . Then, letting $ν ⋆ = \frac{1}{u_{0}} u$ , we obtain that $(x ⋆, ν ⋆)$ is a KKT point. □

Remark 2.20

One of the major difficulties in applying the foregoing result is that we do not know a priori which constraints are active and which constraints are inactive, i.e., the active set is unknown. Therefore, it is necessary to investigate all possible active sets for finding candidate points satisfying the KKT conditions. This is illustrated in the example below.

Example 2.12 (Regular Case). Consider the problem

\begin{matrix} \min_{x \in ℝ^{3}} f (x) : & = \frac{1}{2} (x_{1}^{2} + x_{2}^{2} + x_{3}^{2}) \\ s . t . & g_{1} (x) : = x_{1} + x_{2} + x_{3} + 3 \leq 0 \\ g_{2} (x) : = x_{1} \leq 0 \end{matrix}

Note that every feasible point is regular since $\nabla g_{1} = {[1 1 1]}^{⊤}$ and $\nabla g_{2} = {[1 0 0]}^{⊤}$ , so $x ⋆$ must satisfy the stationarity conditions:

\begin{aligned} x_{1}^{⋆} + ν_{1}^{⋆} + ν_{2}^{⋆} & = 0 \\ x_{2}^{⋆} + ν_{1}^{⋆} & = 0 \\ x_{3}^{⋆} + ν_{1}^{⋆} & = 0 \end{aligned}

Four cases can be distinguished:

The constraints $g_{1}$ and $g_{2}$ are both inactive, i.e. $x_{1}^{⋆} + x_{2}^{⋆} + x_{3}^{⋆} < - 3$ , $x_{1}^{⋆} < 0$ , and $ν_{1}^{⋆} = ν_{2}^{⋆} = 0$ . From the latter together with the dual feasibility conditions, we get $x_{1}^{⋆} = x_{2}^{⋆} = x_{3}^{⋆} = 0$ , hence contradicting the former primal feasibility condition.
The constraint $g_{1}$ is inactive, while $g_{2}$ is active, i.e. $x_{1}^{⋆} + x_{2}^{⋆} + x_{3}^{⋆} < - 3$ , $x_{1}^{⋆} = 0$ , $ν_{1}^{⋆} = 0$ and $ν_{2}^{⋆} \geq 0$ . From the latter, together with the stationarity condition, we get $x_{2}^{⋆} = x_{3}^{⋆} = 0$ , hence contradicting the former once again.
The constraint $g_{1}$ is active, while $g_{2}$ is inactive, i.e. $x_{1}^{⋆} + x_{2}^{⋆} + x_{3}^{⋆} = - 3$ , $x_{1}^{⋆} < 0$ , $ν_{1}^{⋆} \geq 0$ and $ν_{2}^{⋆} = 0$ . Then, the point $(x ⋆, ν ⋆)$ such that $x_{1}^{⋆} = x_{2}^{⋆} = x_{3}^{⋆} = - 1$ , $ν_{1}^{⋆} = 1$ and $ν_{2}^{⋆} = 0$ is a KKT point.
The constraints $g_{1}$ and $g_{2}$ are both active, i.e. $x_{1}^{⋆} + x_{2}^{⋆} + x_{3}^{⋆} = - 3$ , $x_{1}^{⋆} = 0$ , and $ν_{1}^{⋆}, ν_{2}^{⋆} = 0$ . Then, we obtain $x_{2}^{⋆} = x_{3}^{⋆} = - \frac{3}{2}$ , $ν_{1}^{⋆} = \frac{3}{2}$ , and $ν_{2}^{⋆} = - \frac{3}{2}$ , hence contradicting the dual feasibility condition $ν_{2}^{⋆} \geq 0$ .

Overall, there is a unique candidate for a local minimum. Yet, it cannot be concluded as to whether this point is actually a global minimum or even a local minimum.

Remark 2.21: Constraint Qualification

It is very important to note that for a local minimum $x ⋆$ to be a KKT point, an additional condition must be placed on the behaviour of the constraint, i.e., not every local minimum is a KKT point, such a condition is known as a constraint qualification. In Theorem 2.10 it is shown that one possible constraint qualification is that $x ⋆$ be a regular point, which is the well known linear independence constraint qualification (LICQ). A weaker constraint qualification (i.e., implied by LICQ) known as the Mangasarian-Fromovitz constraint qualification (MFCQ) requires that there exits (at least) one direction $d \in 𝒟_{0} (x ⋆)$ , i.e. such that $\nabla g_{i} (x ⋆) ⊤ d < 0$ , for each $i \in 𝒜 (x ⋆)$ . Note, however, that the Lagrange multipliers are guaranteed to be unique if LICQ holds (as stated in Theorem 2.10),while this uniqueness property may be lost under MFCQ.

The following example illustrates the necessity of having a constraint qualification for a KKT point to be a local minimum point of an NLP.

Example 2.13 (Non-Regular Case). Consider the problem

\begin{aligned} \min_{x \in ℝ^{2}} f (x) : & = x_{1}^{2} + {(x_{2} - 1)}^{2} \\ s.t. & g_{1} (x) : = x_{2}^{3} - 2 x_{1} \leq 0 \\ g_{2} (x) : = x_{2}^{3} + 2 x_{1} \leq 0 \end{aligned}

The feasible region is shown in fig. 2.7 below. Note that a minimum point is $x ⋆ = {(0, 0)}^{⊤}$ . The stationary condition relative to variable $x_{2}$ computed at $x ⋆$ reads:

- 2 = 0

It is readily seen that this condition cannot be met at the local minimum point $x ⋆$ . In other words, the KKT conditions are not necessary in this example. This is because no constraint qualification can hold at $x ⋆$ . In particular, $x ⋆$ not being a regular point, LICQ does not hold. Moreover, the set $𝒟_{0} (x ⋆)$ being empty (the direction $d = {(0, - 1)}^{⊤}$ gives $\nabla g_{1} (x ⋆) ⊤ d = \nabla g_{2} (x ⋆) ⊤ d = 0$ while any other direction induces a violation of either one of the constraints), MFCQ does not hold at $x ⋆$ either.

The following theorem provides a sufficient condition under which any KKT point of an inequality constrained NLP problem is guaranteed to be a global minimum of that problem.

Theorem 2.11: KKT Sufficient Conditions

Let $f : ℝ^{n_{x}} \to ℝ$ and $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ , be convex and differentiable functions. Consider the problem to minimize $f (x)$ subject to $g (x) \leq 0$ . If $(x ⋆, ν ⋆)$ is a KKT point, then $x ⋆$ is a global minimum of that problem.

Proof. Consider the function $ℒ (x) : = f (x) + \sum_{i = 1}^{n_{g}} ν_{i}^{⋆} g_{i} (x)$ . Since $f$ and $g_{i}$ , $i = 1, \dots, n_{g}$ are convex functions and $ν_{i} \geq 0$ , $ℒ$ is also convex ^a ^a The sum of convex functions is a convex function, a positive scalar times a convex function is again a convex function . Moreover, the stationarity conditions impose that we have $\nabla ℒ (x ⋆) = 0$ . Hence, by Theorem 2.5, $x ⋆$ is a global minimizer for $ℒ$ on $ℝ^{n_{x}}$ , i.e.

ℒ (x) \geq ℒ (x ⋆) \forall x \in ℝ^{n_{x}}

In particular, for each $x$ such that $g_{i} (x) \leq g_{i} (x ⋆) = 0$ , $i \in 𝒜 (x ⋆)$ , we have

f (x) - f (x ⋆) \geq - \sum_{i \in 𝒜 (x ⋆)} ν_{i}^{⋆} [g_{i} (x) - g_{i} (x ⋆)] \geq 0

Noting that ${x \in ℝ^{n_{x}} : g_{i} (x) \leq 0, i \in 𝒜 (x ⋆)}$ contains the feasible domain ${x \in ℝ^{n_{x}} : g_{i} (x) \leq 0, i = 1, \dots, n_{g}}$ , we therefore showed that $x ⋆$ is a global minimizer for the problem. □

Example 2.14. Consider the problem

\begin{matrix} \min_{x \in ℝ^{3}} f (x) : & = \frac{1}{2} (x_{1}^{2} + x_{2}^{2} + x_{3}^{2}) \\ s . t . & g_{1} (x) : = x_{1} + x_{2} + x_{3} + 3 \leq 0 \\ g_{2} (x) : = x_{1} \leq 0 \end{matrix}

The point $(x ⋆, ν ⋆)$ with $x_{1}^{⋆} = x_{2}^{⋆} = x_{3}^{⋆} = - 1$ , $ν_{1}^{⋆} = 1$ and $ν_{2}^{⋆} = 0$ , being a KKT point, and both the objective function and the feasible set being convex, by Theorem 2.11 , $x ⋆$ is a global minimum.

Both second-order necessary and sufficient conditions for inequality constrained NLP problems will be presented later on in 2.8.

2.7 Problems with equality constraints

In this section, we shall consider nonlinear programming problems with equality constraints of the form:

\begin{matrix} minimize & f (x) \\ subject to & h_{i} (x) = 0 i = 1, \dots, n_{h} \end{matrix}

Based on the material presented in 2.6, it is tempting to convert this problem into an inequality constrained problem, by replacing each equality constraints $h_{i} (x) = 0$ by two inequality constraints $h_{i}^{+} (x) = h_{i} (x) \leq 0$ and $h_{i}^{-} (x) = - h (x) \leq 0$ . Given a feasible point $\bar{x} \in ℝ^{n_{x}}$ , we would have $h_{i}^{+} (\bar{x}) = h_{i}^{-} (\bar{x}) = 0$ and $\nabla h_{i}^{+} (\bar{x}) = - \nabla h_{i}^{-} (\bar{x})$ . Therefore, there could exist no vector $d$ such that $\nabla h_{i}^{+} (\bar{x}) < 0$ and $\nabla h_{i}^{-} (\bar{x}) < 0$ simultaneously, i.e. $𝒟_{0} (\bar{x}) = \emptyset$ . In other words, the geometric conditions developed in the previous section for inequality constrained problems are satisfied by all feasible solutions and, hence, are not informative. A different approach must therefore be used to deal with equality constrained problems. After a number of preliminary results in 2.7.1, we shall describe the method of Lagrange multipliers for equality constrained problems in 2.7.2.

2.7.1 Preliminaries

An equality constraint $h (x) = 0$ defines a set on $ℝ^{n_{x}}$ , which can be seen as a hypersurface.

When considering $n_{h} \geq 1$ equality constraints ${h_{1} (x), \dots, h_{n_{h}} (x)}$ , their intersection forms a (possibly empty) set $S : = {x \in ℝ^{n_{x}} : h_{i} (x) = 0, i = 1, \dots, n_{h}}$ .

Throughout this section, we shall assume that the equality constraints are differentiable; that is, the set ${S : = x \in ℝ^{n_{x}} : h_{i} (x) = 0, i = 1, \dots, n_{h}}$ is said to be a differentiable manifold (or smooth manifold). Associated with a point on a differentiable manifold is the tangent set at that point. To formalize this notion, we start by defining curves on a manifold. A curve $ξ$ on a manifold $S$ is a continuous application $ξ : ℐ \subset ℝ \mapsto S$ , i.e., a family of points $ξ (t) \in S$ continuously parameterized by $t$ in an interval $ℐ \subset ℝ$ . A curve is said to pass through the point $\bar{x}$ if $\bar{x} = ξ (\bar{t})$ for some $\bar{t} \in ℐ$ ; the derivative of a curve at $\bar{t}$ , provided it exists, is defined as $\dot{ξ} (\bar{t}) : = \lim_{h \to 0} \frac{ξ (\bar{t} + h) - ξ (\bar{t})}{h}$ . A curve is said to be differentiable (or smooth) if a derivative exists for each $t \in I$ .

Definition 2.21: Tangent Set

Let $S$ be a (differentiable) manifold in $ℝ^{n_{x}}$ , and let $\bar{x} \in S$ . Consider the collection of all the continuously differentiable curves on $S$ passing through $\bar{x}$ . Then, the collection of all the vectors tangent to these curves at $\bar{x}$ is said to be the tangent set to $S$ at $\bar{x}$ , denoted by $𝒯 (\bar{x})$ .

If the constraints are regular, in the sense of Definition 2.22 below, then $S$ is (locally) of dimension $n_{x} - n_{h}$ , and $𝒯 (\bar{x})$ constitutes a subspace of dimension $n_{x} - n_{h}$ , called the tangent space.

Definition 2.22: Regular Point (for a Set of Equality Constraints)

Let $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ be differentiable functions on $ℝ^{n_{x}}$ and consider the set $S : = {\bar{x} : h_{i} (x) = 0, i = 1, . . ., n_{h}}$ . A point $x \in S$ is said to be a regular point if the gradient vectors $\nabla h_{i} (x)$ , $i = 1, \dots, n_{h}$ are linearly independent, i.e.,

rank (\nabla h_{1} (\bar{x}) \nabla h_{2} (\bar{x}) \dots \nabla h_{n_{h}} (\bar{x})) = n_{h}

Lemma 2.3: Algebraic Characterization of a Tangent Space

Let $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ be differentiable functions on $ℝ^{n_{x}}$ and consider the set $S : = {\bar{x} : h_{i} (x) = 0, i = 1, . . ., n_{h}}$ . At a regular point $x \in S$ , the tangent space is such that

𝒯 (\bar{x}) = {d : \nabla h {(\bar{x})}^{⊤} d = 0}

Proof. ... □

Recall that $\nabla h ⊤ = \frac{\partial h}{\partial x} \in ℝ^{n_{h} \times n_{x}}$ is the Jacobian of the constraints. Therefore, the tangent directions $d$ must be orthogonal to the gradient of each constraint $\nabla h_{i}$ at $\bar{x}$ . Note also that $\frac{\partial h}{\partial x} d = 0$ defines the kernel of the Jacobian at $\bar{x}$ that is the set $K (\frac{\partial h}{\partial x}) : = {d : \frac{\partial h}{\partial x} d = 0}$ . Since we have assumed that $\bar{x}$ is a regular point, the Jacobian is full rank and for the rank-nullity theorem we have that the dimension of the kernel is $dim (K (\frac{\partial h}{\partial x})) = n_{x} - n_{h}$ that is the dimension of the tangent space.

2.7.2 The Method of Lagrange Multipliers

The idea behind the method of Lagrange multipliers for solving equality constrained NLP problems of the form

\begin{matrix} minimize & f (x) \\ subject to & h_{i} (x) = 0 i = 1, \dots, n_{h} \end{matrix}

is to restrict the search of a minimum on the manifold $S : = {x \in ℝ^{n_{x}} : h_{i} (x) = 0, i = 1, . . ., n_{h}}$ . In other words, we derive optimality conditions by considering the value of the objective function along curves on the manifold S passing through the optimal point.

The following Theorem shows that the tangent space $𝒯 (x ⋆)$ at a regular (local) minimum point $x ⋆$ is orthogonal to the gradient of the objective function at $x ⋆$ . This important fact is illustrated in Fig. ??. in the case of a single equality constraint.

Theorem 2.12: Geometric Necessary Condition for a Local Minimum

Let $f : ℝ^{n_{x}} \to ℝ$ and $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ be continuously differentiable functions on $ℝ^{n_{x}}$ . Suppose that $x ⋆$ is a local minimum point of the problem to minimize $f (x)$ subject to the constraints $h (x) = 0$ . Then, $\nabla f (x ⋆)$ is orthogonal to the tangent space $𝒯 (x ⋆)$ , that is:

ℱ_{0} (x ⋆) \cap 𝒯 (x ⋆) = \emptyset

Proof. By contradiction, assume that there exists a $d \in 𝒯 (x ⋆)$ such that $\nabla f {(x}^{⋆}) ⊤ d \neq 0$ . Let $ξ : ℐ = [- a, a] \to S, a > 0$ be any smooth curve passing through $x ⋆$ with $ξ (0) = x ⋆$ and $\dot{ξ} (0) = d$ . Let also $ϕ$ be the function defined as $ϕ (t) : = f (ξ (t)), \forall t \in I$ . Since $x ⋆$ is a local minimum of f on $S : = {x : h_{i} (x) = 0, i = 1, . . ., n_{h}}$ , by Definition 2.11, we have

\exists η > 0 such that φ (t) = f (ξ (t)) \geq f (x ⋆) = φ (0) \forall t \in ℬ_{η} (0) \cap ℐ

It follows that $t^{⋆} = 0$ is an unconstrained (local) minimum point for $ϕ$ , and

0 = \nabla φ (0) = \nabla f (x ⋆) ⊤ \dot{ξ} (0) = \nabla f (x ⋆) ⊤ d

We thus get a contradiction with the fact that $\nabla f {(x}^{⋆}) ⊤ d \neq 0$ □

Next, we take advantage of the forgoing geometric characterization, and derive first-order necessary conditions for equality constrained NLP problems.

Theorem 2.13: First-Order Necessary Conditions

Let $f : ℝ^{n_{x}} \to ℝ$ and $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ be continuously differentiable functions on $ℝ^{n_{x}}$ . Consider the problem to minimize $f (x)$ subject to the constraints $h (x) = 0$ . If $x ⋆$ is a local minimum and is a regular point of the constraints, then there exists a unique vector $λ ⋆ \in ℝ^{n_{h}}$ such that:

\nabla f (x ⋆) + \nabla h (x ⋆) λ ⋆ = \nabla f (x ⋆) + \sum_{i = 1}^{n_{h}} \nabla h_{i} (x ⋆) λ_{i}^{⋆} = 0

Proof. Since $x ⋆$ is a local minimum of $f$ on $S : = {x \in ℝ^{n_{x}} : h (x) = 0}$ , by Theorem 2.12, we have $ℱ_{0} (x ⋆) \cap 𝒯 (x ⋆) = \emptyset$ , i.e. the system

\nabla f (x ⋆) ⊤ d < 0 \nabla h (x ⋆) ⊤ d = 0

is inconsistent. Consider the following two sets:

\begin{matrix} C_{1} : = {(z_{1}, z 2) \in ℝ^{n_{h} + 1} : z_{1} = \nabla f (x ⋆) ⊤ d, z 2 = \nabla h (x ⋆) ⊤ d} \\ C_{2} : = {(z_{1}, z 2) \in ℝ^{n_{h} + 1} : z_{1} < 0, z 2 = 0 .} \end{matrix}

Clearly, $C_{1}$ and $C_{2}$ are convex and $C_{1} \cap C_{2} = \emptyset$ . Then, by the separation Theorem 2.14, there exists a nonzero vector $(μ, λ) \in ℝ^{n_{h} + 1}$ such that

μ \nabla f (x ⋆) ⊤ d + λ ⊤ [\nabla h (x ⋆) ⊤ d] \geq μ z_{1} + λ ⊤ z 2 \forall d \in ℝ^{n_{x}}, \forall (z_{1}, z 2) \in cl (C_{2})

Letting $z 2 = 0$ and since $z_{1}$ can be made an arbitrarily large negative number, it follows that $μ \geq 0$ ^a ^a Note that if $μ < 0$ then since $z_{1} < 0$ the lower bound on the left-hand side would be positive and arbitrarily large. Also, letting $(z_{1}, z 2) = (0, 0)$ ^b we must have $[μ \nabla f (x ⋆) + \nabla h (x ⋆) λ] ⊤ d \geq 0$ , for each $d \in ℝ^{n_{x}}$ . In particular, letting $d = - [μ \nabla f (x ⋆) + \nabla h (x ⋆) λ]$ , it follows that $- ∥ μ \nabla f (x ⋆) + \nabla h (x ⋆) λ ∥ 2 \geq 0$ , and thus,

μ \nabla f (x ⋆) + \nabla h (x ⋆) λ = 0 with (μ, λ) \neq (0, 0)

Finally, note that $μ > 0$ , for otherwise the above equation would contradict the assumption of linear independence of $\nabla h_{i} (x ⋆), i = 1, \dots, n_{h}$ . The result follows by letting $λ ⋆ : = \frac{1}{μ} λ$ , and noting that the linear independence assumption implies the uniqueness of these Lagrangian multipliers. □

Theorem 2.14: Separation of Two Convex Sets

Let $C_{1}$ and $C_{2}$ be two nonempty, convex set in $ℝ^{n}$ and suppose that $C_{1} \cap C_{2} = \emptyset$ . Then, there exists a hyperplane that separates $C_{1}$ and $C_{2}$ ; that is, there exists a nonzero vector $p \in ℝ^{n}$ such that

p ⊤ x 1 \geq p ⊤ x 2 \forall x 1 \in cl (C_{1}) \forall x 2 \in cl (C_{2})

Remark 2.22: Obtaining Candidate Solution Points

The first-order necessary conditions

\nabla f (x ⋆) + \nabla h (x ⋆) λ ⋆ = 0

together with the constraints

h (x ⋆) = 0

give a total of $n_{x} + n_{h}$ (typically nonlinear) equations in the variables $(x ⋆, λ ⋆)$ . Hence, these conditions are complete in the sense that they determine, at least locally, a unique solution. However, as in the unconstrained case, a solution to the first-order necessary conditions need not be a (local) minimum of the original problem; it could very well correspond to a (local) maximum point, or some kind of saddle point. These consideration are illustrated in Example 2.15 below.

Remark 2.23: Regularity-Type Assumption

It is important to note that for a local minimum to satisfy the foregoing first-order conditions and, in particular, for a unique Lagrange multiplier vector to exist, it is necessary that the equality constraint satisfy a regularity condition. In other word, the first-order conditions may not hold at a local minimum point that is non-regular. An illustration of these considerations is provided in Example 2.16.

There exists a number of similarities with the constraint qualification needed for a local minimizer of an inequality constrained NLP problem to be KKT point; in particular, the condition that the minimum point be a regular point for the constraints corresponds to LICQ (see Remark 2.21).

Remark 2.24: Lagrangian

It is convenient to introduce the Lagrangian $ℒ : ℝ^{n_{x}} \times ℝ^{n_{h}} \to ℝ$ associated with the constrained problem, by adjoining the cost and constraint functions as:

ℒ (x, λ) : = f (x) + λ ⊤ h (x)

Thus, if $x ⋆$ is a local minimum which is regular, the first-order necessary conditions are written as

\begin{matrix} \nabla_{x} ℒ (x ⋆, λ ⋆) = 0 \\ \nabla_{λ} ℒ (x ⋆, λ ⋆) = 0 \end{matrix}

the latter equations being simply a restatement of the constraints. Note that the solution of the original problem typically corresponds to a saddle point of the Lagrangian function.

Example 2.15 (Regular Case). Consider the problem

\begin{aligned} \min_{x \in ℝ^{2}} f (x) & : = x_{1} + x_{2} \\ s.t. g (x) & : = x_{1}^{2} + x_{2}^{2} - 2 = 0 \end{aligned}

Observe first that every feasible point is a regular point for the equality constraint (the point (0,0) being infeasible). Therefore, every local minimum is a stationary point of the Lagrangian function by Theorem 2.13 .

The gradient vectors $\nabla f (x)$ and $\nabla h (x)$ are given by

\nabla f (x) = {(\begin{matrix} 1 & 1 \end{matrix})}^{⊤} and \nabla h (x) = {(\begin{matrix} 2 x_{1} & 2 x_{2} \end{matrix})}^{⊤}

so that the first-order necessary conditions read

\begin{aligned} 2 λ x_{1} & = - 1 \\ 2 λ x_{2} & = - 1 \\ x_{1}^{2} + x_{2}^{2} & = 2 \end{aligned}

These three equations can be solved for the three unknowns $x_{1}$ , $x_{2}$ and $λ$ . Two candidate local minimum points are obtained: (i) $x_{1}^{⋆} = x_{2}^{⋆} = - 1, λ^{⋆} = 1$ , and (ii) $x_{1}^{⋆} = x_{2}^{⋆} = 1, λ^{⋆} = - 1$ . These results are illustrated on Fig. 2.8 . It can be seen that only the former actually corresponds to a local minimum point, while the latter gives a local maximum point.

Example 2.16 (Non-Regular Case). Consider the problem

\begin{matrix} \min_{x \in ℝ^{2}} & f (x) : = - x_{1} \\ s . t . & h_{1} (x) : = {(1 - x_{1})}^{3} + x_{2} = 0 \\ h_{2} (x) : = {(1 - x_{1})}^{3} - x_{2} = 0 \end{matrix}

As shown by Fig. 1.13., this problem has only one feasible point, namely, $x ⋆ = {(1 0)}^{⊤}$ ; that is, $x ⋆$ is also the unique global minimum of the problem. However, at this point, we have

\nabla f (x ⋆) = (\begin{matrix} - 1 \\ 0 \end{matrix}), \nabla h_{1} (x ⋆) = (\begin{matrix} 0 \\ 1 \end{matrix}) and \nabla h_{2} (x ⋆) = (\begin{matrix} 0 \\ - 1 \end{matrix})

hence the first-order conditions

λ_{1} (\begin{matrix} 0 \\ 1 \end{matrix}) + λ_{2} (\begin{matrix} 0 \\ - 1 \end{matrix}) = (\begin{matrix} 1 \\ 0 \end{matrix})

cannot be satisfied. This illustrates the fact that a minimum point may not be a stationary point for the Lagrangian if that point is non-regular.

The following theorem provides second-order necessary conditions for a point to be a local minimum of a NLP problem with equality constraints.

Theorem 2.15: Second-Order Necessary Conditions

Let $f : ℝ^{n_{x}} \to ℝ$ and $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ be twice continuously differentiable functions on $ℝ^{n_{x}}$ . Consider the problem to minimize $f (x)$ subject to the constraints $h (x) = 0$ . If $x ⋆$ is a local minimum and is a regular point of the constraints, then there exists a unique vector $λ ⋆ \in ℝ^{n_{h}}$ such that

\nabla f (x ⋆) + \nabla h (x ⋆) λ ⋆ = 0

and

d ⊤ (\nabla^{2} f (x ⋆) + \sum_{i = 1}^{n_{h}} \nabla^{2} h_{i} (x ⋆) λ_{i}^{⋆}) d \geq 0 \forall d such that \nabla h (x ⋆) ⊤ d = 0

Proof. Note first that $\nabla f (x ⋆) + \nabla h (x ⋆) λ ⋆ = 0$ directly follows from Theorem 2.13. Let $d$ be an arbitrary direction in $𝒯 (x ⋆)$ ; that is, $\nabla h {(x}^{⋆}) ⊤ d = 0$ since $x ⋆$ is a regular point (see Lemma 2.3). Consider any twice differentiable curve $ξ : ℐ = [- a, a] \to S$ , $a > 0$ passing through $x ⋆$ with $ξ (0) = x ⋆$ and $\dot{ξ} (0) = d$ . Let $ϕ$ be the function defined as $ϕ (t) : = f (ξ (t)) \forall t \in ℐ$ . Since $x ⋆$ is a local minimum of $f$ on $S : = {x \in ℝ^{n_{x}} : h (x) = 0}$ , $t^{⋆} = 0$ is an unconstrained (local) minimum point for $ϕ$ . By Theorem 2.4, it follows that

0 \leq \nabla^{2} φ (0)

^a ^a Note that ϕ is a scalar function, as a consequence ∇ ⁡ 2 φ = d 2 ϕ d t 2 = ξ ˙ ( 0 ) ⊤ ⁡ ∇ ⁡ 2 f ( x ⋆ ) ξ ˙ ( 0 ) + ∇ ⁡ f ( x ⋆ ) ⊤ ⁡ ξ ¨ ( 0 )

Furthermore, differentiating the relation $h {(ξ (t))}^{⊤} λ ⋆ = 0$ twice, we obtain

\dot{ξ} {(0)}^{⊤} (\sum_{i = 1}^{n_{h}} \nabla^{2} h_{i} (x ⋆) λ_{i}^{⋆}) \dot{ξ} (0) + (\nabla h (x ⋆) λ ⋆) ⊤ \ddot{ξ} (0) = 0

Adding the last two equations yields

d ⊤ (\nabla^{2} f (x ⋆) + \sum_{i = 1}^{n_{h}} \nabla^{2} h_{i} (x ⋆) λ_{i}^{⋆}) d \geq 0

and this condition must hold for every $d$ such that $\nabla h (x ⋆)^{⊤} d = 0$ □

It is useful to shed more light on the derivation of the previous proof. Thus, we carry out the derivative of $ψ (t) = h {(ξ (t))}^{⊤} λ ⋆$ explicitly. Note that since the curve $ξ$ is in $S$ we have $h_{i} (ξ (t)) = 0 \forall i = 1, \dots, n_{h} \forall t \in ℐ$ . Hence $ψ (t)$ is identically zero for all $t$ (i.e. $ψ (t) \equiv 0 \forall t \in ℐ$ ). We can expand $ψ$ as $ψ (t) = \sum_{i = 1}^{n_{h}} h_{i} (ξ (t)) λ_{i}^{⋆}$ . The first total derivative is:

\dot{ψ} = \sum_{i = 1}^{n_{h}} \sum_{j = 1}^{n_{x}} \frac{\partial h_{i}}{\partial x_{j}} {\dot{ξ}}_{j} λ_{i}^{⋆} = [\nabla h λ ⋆] ⊤ \dot{ξ} = λ ⊤ \frac{\partial h}{\partial x} \dot{ξ} = 0

We can now differentiate again with respect to time to derive the second total derivative:

\begin{aligned} \ddot{ψ} = \sum_{i = 1}^{n_{h}} \sum_{j = 1}^{n_{x}} \sum_{k = 1}^{n_{x}} \frac{\partial^{2} h_{i}}{\partial x_{k} \partial x_{j}} {\dot{ξ}}_{k} {\dot{ξ}}_{j} λ_{i}^{⋆} + \sum_{i = 1}^{n_{h}} \sum_{j = 1}^{n_{x}} \frac{\partial h_{i}}{\partial x_{j}} {\ddot{ξ}}_{j} λ_{i}^{⋆} = \\ = \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} {\dot{ξ}}^{⊤} \nabla^{2} h_{i} \dot{ξ} + \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} \nabla h_{i}^{⊤} \ddot{ξ} = 0 \end{aligned}

Remember that in $ξ (0) = x ⋆$ and $\dot{ξ} (0) = d$ and the first-order necessary condition hold $\nabla f (x ⋆) = - \nabla h (x ⋆) λ ⋆ = - \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} \nabla h_{i} (x ⋆)$ . Therefore, the second derivative of $ψ$ in $0$ gives:

\ddot{ψ} (0) = \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} d ⊤ \nabla^{2} h_{i} (x ⋆) d - \nabla f {(x}^{⋆}) ⊤ \ddot{ξ} (0) = 0

Hence:

\nabla f {(x}^{⋆}) ⊤ \ddot{ξ} (0) = \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} d ⊤ \nabla^{2} h_{i} (x ⋆) d

Using the last equation, it is easier to follow the proof of Theorem 2.15

Remark 2.25: Eigenvalues in Tangent Space

In the foregoing Theorem, it is shown that the matrix $\nabla x x 2 ℒ (x ⋆, λ ⋆)$ restricted to the subspace $𝒯 (x ⋆)$ plays a key role. Geometrically, the restriction of $\nabla x x 2 ℒ (x ⋆, λ ⋆)$ to $𝒯 (x ⋆)$ corresponds to the projection
$𝒫_{𝒯 (x ⋆)} [\nabla x x 2 ℒ (x ⋆, λ ⋆)]$ .

A vector $y \in 𝒯 (x ⋆)$ is said to be an eigenvector of $𝒫_{𝒯 (x ⋆)} [\nabla x x 2 ℒ (x ⋆, λ ⋆)]$ if there is a real number $μ$ such that:

𝒫_{𝒯 (x ⋆)} [\nabla x x 2 ℒ (x ⋆, λ ⋆)] y = μ y

the corresponding $μ$ is said to be an eigenvalue of $𝒫_{𝒯 (x ⋆)} [\nabla x x 2 ℒ (x ⋆, λ ⋆)]$ ^a ^a These definitions coincide with the usual definitions of eigenvector and eigenvalue for real matrices..

Now, to obtain a matrix representation for $𝒫_{𝒯 (x ⋆)} [\nabla x x 2 ℒ (x ⋆, λ ⋆)]$ , it is necessary to introduce a basis of the subspace $𝒯 (x ⋆)$ , say
$E = [e 1, . . ., e n_{x} - n_{h}]$ ^b ^b Note that $E \in ℝ^{n_{x} \times (n_{x} - n_{h})}$ . Then, the eigenvalues of $𝒫_{𝒯 (x ⋆)} [\nabla x x 2 ℒ (x ⋆, λ ⋆)]$ are the same as those of the $n_{x} - n_{h} \times n_{x} - n_{h}$ matrix $E^{⊤} \nabla x x 2 ℒ (x ⋆, λ ⋆) E$ ; in particular, they are independent of the particular choice of basis $E$ .

Example 2.17 (Regular Case Continued). Consider again the Example 2.15 . Two candidate local minimum points, (i) $x_{1}^{⋆} = x_{2}^{⋆} = - 1, λ^{⋆} = 1$ and (ii) $x_{1}^{⋆} = x_{2}^{⋆} = 1, λ^{⋆} = - 1$ , were obtained on application of the first-order necessary conditions. The Hessian matrix of the Lagrangian function is given by

\nabla_{x x}^{2} ℒ (x, λ) = \nabla^{2} f (x) + λ \nabla^{2} h (x) = λ (\begin{matrix} 2 & 0 \\ 0 & 2 \end{matrix})

and a basis of the tangent subspace at a point $x \in 𝒯 (x)$ is

E (x) : = (\begin{matrix} - x_{2} \\ x_{1} \end{matrix})

Therefore,

E^{⊤} \nabla_{x x}^{2} ℒ (x, λ) E = 2 λ (x_{1}^{2} + x_{2}^{2})

Since all the feasible $x_{1}$ and $x_{2}$ must satisfy the constraint $h = x_{1}^{2} + x_{2}^{2} - 2 = 0$ we have:

E^{⊤} \nabla_{x x}^{2} ℒ (x, λ) E = 4 λ

In particular, for the candidate solution point (i), we have

E^{⊤} \nabla_{x x}^{2} ℒ (x ⋆, λ^{⋆}) E = 4 > 0

hence satisfying the second-order necessary conditions (in fact, this point also satisfies the second-order sufficient conditions of optimality discussed hereafter). On the other hand, for the candidate solution point (ii), we get

E^{⊤} \nabla_{x x}^{2} ℒ (x ⋆, λ^{⋆}) E = - 4 < 0

which does not satisfy the second-order necessary requirement, so this point cannot be a local minimum.

The conditions given in Theorems 2.13 and 2.15 are necessary conditions that must hold at each local minimum point. Yet, a point satisfying these conditions may not be a local minimum. The following theorem provides sufficient conditions for a stationary point of the Lagrangian function to be a (local) minimum, provided that the Hessian matrix of the Lagrangian function is locally convex along directions in the tangent space of the constraints.

Theorem 2.16: Second-Order Sufficient Conditions

\begin{matrix} \nabla_{x} ℒ (x ⋆, λ^{⋆}) = \nabla f (x ⋆) + \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} \nabla h_{i} (x ⋆) = 0 \\ \nabla_{λ} ℒ (x ⋆, λ^{⋆}) = h (x ⋆) = 0 \end{matrix}

and

y ⊤ \nabla_{x x}^{2} ℒ (x ⋆, λ ⋆) y > 0 \forall y \neq 0 such that \nabla h (x ⋆) ⊤ y = 0

where $ℒ (x, λ) = f (x) + λ ⊤ h (x)$ , then $x ⋆$ is a strict local minimum.

Proof. Consider the augmented Lagrangian function

\bar{ℒ} (x, λ) = f (x) + λ ⊤ h (x) + \frac{c}{2} ∥ h (x) ∥^{2}

where $c$ is a scalar. We have

\begin{matrix} \nabla_{x} \bar{ℒ} (x, λ) & = \nabla_{x} ℒ (x, \bar{λ}) \\ \nabla_{x x}^{2} \bar{ℒ} (x, λ) & = \nabla_{x x}^{2} ℒ (x, \bar{λ}) + c \nabla h {(x)}^{⊤} \nabla h (x) \end{matrix}

where $\bar{λ} = λ + c h (x)$ . Since $(x ⋆, λ ⋆)$ satisfy the sufficient conditions and by Lemma 2.4, we obtain

\nabla_{x} \bar{ℒ} (x ⋆, λ ⋆) = 0 and \nabla_{x x}^{2} \bar{ℒ} (x ⋆, λ ⋆) ≽ 0

for sufficiently large $c$ . $\bar{ℒ}$ being definite positive at $(x ⋆, λ ⋆)$ ,

\begin{aligned} \exists ϱ > 0, δ > 0 such that \bar{ℒ} (x, λ ⋆) \geq \bar{ℒ} (x ⋆, λ ⋆) + \frac{ϱ}{2} ∥ x - x ⋆ ∥ 2 \\ for ∥ x - x ⋆ ∥ < δ \end{aligned}

Finally, since $\bar{ℒ} (x ⋆, λ ⋆) = f (x)$ when $h (x) = 0$ , we get

f (x) \geq f (x ⋆) + \frac{ϱ}{2} ∥ x - x ⋆ ∥ 2 if h (x) = 0, ∥ x - x ⋆ ∥ < δ

i.e., $x ⋆$ is a strict local minimum. □

Example 2.18. Consider the problem

\begin{matrix} \min_{x \in ℝ^{3}} f (x) : = - x_{1} x_{2} - x_{1} x_{3} - x_{2} x_{3} \\ s . t . h (x) : = x_{1} + x_{2} + x_{3} - 3 = 0 \end{matrix}

The first-order conditions for this problem are

\begin{matrix} - (x_{2} + x_{3}) + λ = 0 \\ - (x_{1} + x_{3}) + λ = 0 \\ - (x_{1} + x_{2}) + λ = 0 \end{matrix}

together with the equality constraint. It is easily checked that the point $x_{1}^{⋆} = x_{2}^{⋆} = x_{3}^{⋆} = 1$ , $λ^{⋆} = 2$ satisfies these conditions. Moreover,

\nabla_{x x}^{2} ℒ (x ⋆, λ ⋆) = \nabla^{2} f (x ⋆) = (\begin{matrix} 0 & - 1 & - 1 \\ - 1 & 0 & - 1 \\ - 1 & - 1 & 0 \end{matrix})

and a basis of the tangent space to the constraint $h (x) = 0$ at $x ⋆$ is

E : = (\begin{matrix} 0 & 2 \\ 1 & - 1 \\ - 1 & - 1 \end{matrix})

We thus obtain

E ⊤ \nabla_{x x}^{2} ℒ (x ⋆, λ^{⋆}) E = (\begin{matrix} 2 & 0 \\ 0 & 2 \end{matrix})

which is clearly a definite positive matrix. Hence, $x ⋆$ is a strict local minimum of (1.18). (Interestingly enough, the Hessian matrix of the objective function itself is indefinite at $x ⋆$ in this case.)

We close this section by providing insight into the Lagrange multipliers.

Remark 2.26: A first Interpretation of the Lagrange Multipliers

The concept of Lagrange multipliers allows to adjoin the constraints to the objective function. That is, one can view constrained optimization as a search for a vector $x ⋆$ at which the gradient of the objective function is a linear combination of the gradients of constraints.

Another insightful interpretation of the Lagrange multipliers is as follows. Consider the set of perturbed problems $v^{⋆} (y) : = min {f (x) : h (x) = y}$ . Suppose that there is a unique regular solution point for each $y$ , and let $ξ ⋆ (y) : = arg min {f (x) : h (x) = y}$ denote the evolution of the optimal solution point as a function of the perturbation parameter $y$ . Clearly,

v (0) = f (x ⋆) and ξ ⋆ (0) = x ⋆

Moreover, since $h (ξ ⋆ (y)) = y$ for each $y$ , we have:

\nabla_{y} h (ξ ⋆ (y)) = \frac{d h}{d y} = 1 = \sum_{i = 1}^{n_{x}} \frac{\partial h}{\partial x_{j}} \frac{d ξ_{j}^{⋆}}{d y} = \nabla h^{⊤} \frac{d ξ}{⋆}

Denoting by $λ^{⋆}$ the Lagrange multiplier associated to the constraint $h (x) = 0$ in the original problem, we have

\frac{d v^{⋆}}{d y} |_{y = 0} = \nabla f {(x}^{⋆}) ⊤ \frac{d ξ}{⋆} = - λ^{⋆} \nabla h {(x}^{⋆}) ⊤ \frac{d ξ}{⋆} = - λ^{⋆}

Therefore, the Lagrange multipliers $λ^{⋆}$ can be interpreted as the sensitivity of the objective function $f$ with respect to the constraint $h$ . Said differently, $λ^{⋆}$ indicates how much the optimal cost would change, if the constraints were perturbed. This interpretation extends straightforwardly to NLP problems having inequality constraints. The Lagrange multipliers of an active constraints $g (x) \leq 0$ , say $ν^{⋆} > 0$ , can be interpreted as the sensitivity of $f (x ⋆)$ with respect to a change in that constraints, as $g (x) \leq y$ ; in this case, the positivity of the Lagrange multipliers follows from the fact that by increasing $y$ , the feasible region is relaxed, hence the optimal cost cannot increase. Regarding inactive constraints, the sensitivity interpretation also explains why the Lagrange multipliers are zero, as any infinitesimal change in the value of these constraints leaves the optimal cost value unchanged.

Remark 2.27: A second Interpretation of the Lagrange Multipliers: Mechanical Analogy

It is possible to introduce an interesting physical interpretation of the KKT condition and its Lagrange multipliers. Consider a ball on a hilly terrain where the elevation of the terrain is described by the function $h = f (x)$ where $x = {[\begin{matrix} x_{1} & x_{2} \end{matrix}]}^{⊤}$ are the $x$ and $y$ coordinates. The gravitational potential is $V = m g f (x)$ , we use normalized units so that $V = f (x)$ . In the unconstrained case the equilibrium point is a stationary point for the potential that is $\nabla V (x ⋆) = 0$ . The generalized force, due to gravity, acting on the ball in a generic point $x$ is $F = - \nabla V (x) = - \nabla f (x)$ . Note that the motion along the vertical axis is completely determined by $x$ (i.e. the system has two degrees of freedom). An equality constraint can be seen as a rail on which the ball is forced to slide. The rail is modeled as a curve in the form $h (x) = 0$ . Such a constraint exert a generalized reaction force at each point normal to the curve, hence in the direction of $\nabla h (x)$ . The magnitude of this force depends on how much the gravitational force is pushing against the constraint. The static equilibrium condition is $- \nabla f (x ⋆) = λ^{⋆} \nabla h (x ⋆)$ . On the other hand, inequality constraints can be seen as barriers, the ball is forced to lie inside these barriers. Consider the inequality constraint $g (x) \leq 0$ . When the constraint is inactive, it is not exerting any force on the ball and hence it is not affecting the equilibrium equation. When the constraint is active, the direction of the force exerted is $\nabla g (x)$ .Compared to rails, barriers are unilateral constraints and can only exert forces in one sense. This is the physical interpretation of the nonnegativity of its associated Lagrangian multipliers $ν \geq 0$ . The equilibrium equation in presence of both equality and inequality constraints become the first KKT condition. This interpretation is illustrated in Figure 2.9. Gravitational and reaction forces are $F = \nabla f (x ⋆)$ , $R g = - ν^{⋆} \nabla g (x ⋆)$ , $R h = - λ^{⋆} \nabla h (x ⋆)$ . The force balance is therefore the KKT stationary condition with reversed sign : $F + R g + R h = 0$

2.8 General NLP problems

In this section, we shall consider general, nonlinear programming problems with both equality and inequality constraints,

\begin{matrix} minimize & f (x) \\ subject to & g_{i} (x) \leq 0, i = 1, \dots, n_{g} \\ h_{i} (x) = 0, i = 1, \dots, n_{h} \end{matrix}

Before stating necessary and sufficient conditions for such problems, we shall start by revisiting the KKT conditions for inequality constrained problems, based on the method of Lagrange multipliers described in 2.7.

2.8.1 KKT Conditions for Inequality Constrained NLP Problems Revisited

Consider the problem to minimize a function $f (x) for x \in S : = {x \in ℝ^{n_{x}} : g_{i} (x) = 0, i = 1, . . ., n_{g}}$ and suppose that $x ⋆$ is a local minimum point. Clearly, $x ⋆$ is also a local minimum of the inequality constrained problem where the inactive constraints $g_{i} (x) \leq 0, i \notin 𝒜 (x ⋆)$ , have been discarded. Thus, in effect, inactive constraints at $x ⋆$ don’t matter; they can be ignored in the statement of optimality conditions.

On the other hand, active inequality constraints can be treated to a large extent as equality constraints at a local minimum point. In particular, it should be clear that $x ⋆$ is also a local minimum to the equality constrained problem

\begin{matrix} minimize & f (x) \\ subject to & g_{i} (x) = 0, i \in 𝒜 (x ⋆) \end{matrix}

That is, it follows from Theorem 2.13 that, if $x ⋆$ is a regular point, there exists a unique Lagrange multiplier vector $ν ⋆ \in ℝ^{n_{g}}$ such that

\nabla f (x ⋆) + \sum_{i \in 𝒜 (x ⋆)} ν_{i}^{⋆} \nabla g_{i} (x ⋆) = 0

Assigning zero Lagrange multipliers to the inactive constraints, we obtain

\begin{aligned} \nabla f (x ⋆) + \nabla g (x ⋆) ν ⋆ & = 0 \\ ν_{i}^{⋆} & = 0, \forall i \notin 𝒜 (x ⋆) \end{aligned}

This latter condition can be rewritten by means of the following equations:

ν_{i}^{⋆} g_{i} (x ⋆) = 0 \forall i = 1, \dots, n_{g}

The argument showing that $ν ⋆ \geq 0$ is a little more elaborate. By contradiction, assume that $ν_{l} < 0$ for some $l \in 𝒜 (x ⋆)$ . Let $A \in ℝ^{(n_{g} + 1) \times n_{x})}$ be the matrix whose rows are $\nabla f (x ⋆)$ and $\nabla g_{i} (x ⋆), i = 1, . . ., n g$ . Since $x ⋆$ is a regular point, the Lagrange multiplier vector $ν ⋆$ is unique. Therefore, the condition

A^{⊤} y = 0

can only be satisfied by $y ⋆ : = η {(1 ν}^{⋆}) ⊤$ with $η \in ℝ$ . Because $ν_{l} < 0$ , we know by Theorem 2.1 that there exists a direction $\bar{d} \in ℝ^{n_{x}}$ such that $A \bar{d} < 0$ . In other words,

\bar{d} \in ℱ_{0} (x ⋆) \cap 𝒟_{0} (x ⋆) \neq \emptyset

which contradicts the hypothesis that $x ⋆$ is a local minimizer of $f$ on $S$ (see Remark 2.18).

Overall, these results thus constitute the KKT optimality conditions as stated in Theorem 2.10. But although the foregoing development is straightforward, it is somewhat limited by the regularity-type assumption at the optimal solution. Obtaining more general constraint qualifications (see Remark 2.21) requires that the KKT conditions be derived using an alternative approach, e.g., the approach described earlier in Section 2.6. Still, the conversion to equality constrained problem proves useful in many situations, e.g., for deriving second-order sufficiency conditions for inequality constrained NLP problems.

2.8.2 Optimality Conditions for General NLP Problems

We are now ready to generalize the necessary and sufficient conditions given in Theorems 2.10, 2.13, 2.15 and 2.16 to general NLP problems.

Theorem 2.17: First- and Second-Order Necessary Conditions

Let $f : ℝ^{n_{x}} \to ℝ$ , $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ and $h_{i} : ℝ^{n_{x}} \to ℝ$ , be twice continuously differentiable functions on $ℝ^{n_{x}}$ . Consider the problem of minimizing $f (x)$ subject to the constraints $h (x) = 0$ and $g (x) = 0$ . If $x ⋆$ is a local minimum of the optimization problem and is a regular point of the constraints, then there exist unique vectors $ν ⋆$ and $λ ⋆$ such that

\begin{aligned} \nabla f (x ⋆) + \nabla h (x ⋆) λ ⋆ + \nabla g (x ⋆) ν ⋆ & = 0 \\ g (x ⋆) & \leq 0 \\ h (x ⋆) & = 0 \\ g {(x}^{⋆}) ⊤ ν ⋆ & = 0 \end{aligned}

and

y ⊤ (\nabla^{2} f (x ⋆) + \sum_{i = 1}^{n_{g}} ν_{i}^{⋆} \nabla^{2} g_{i} (x ⋆) + \sum_{i = 1}^{n_{h}} λ_{i}^{⋆} \nabla^{2} h_{i} (x ⋆)) y \geq 0

for all $y$ such that $\nabla g_{i} {(x}^{⋆}) ⊤ y = 0, i = 1 \in 𝒜 (x ⋆)$ and $\nabla h {(x}^{⋆}) ⊤ y = 0$

Theorem 2.18: Second-Order Sufficient Conditions

Let $f : ℝ^{n_{x}} \to ℝ$ , $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ and $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ , be twice continuously differentiable functions on $ℝ^{n_{x}}$ . Consider the problem of minimizing $f (x)$ subject to the constraints $h (x) = 0$ and $g (x) = 0$ . If there exists $x ⋆$ , $ν ⋆$ and $λ ⋆$ satisfying the KKT conditions in Theorem 2.17, and

y ⊤ \nabla_{x x}^{2} ℒ (x ⋆, ν ⋆, λ ⋆) y > 0

for all $y \neq 0$ such that

\begin{matrix} \nabla g_{i} {(x}^{⋆}) ⊤ y & = 0 & i \in 𝒜 (x ⋆) with ν_{i}^{⋆} > 0 \\ \nabla g_{i} {(x}^{⋆}) ⊤ y & \leq 0 & i \in 𝒜 (x ⋆) with ν_{i}^{⋆} = 0 \\ \nabla h {(x}^{⋆}) ⊤ y & = 0 \end{matrix}

where $ℒ = f (x) + h {(x)}^{⊤} λ + g {(x)}^{⊤} ν$ , then $x ⋆$ is a strict local minimum.

Likewise, the KKT sufficient conditions given in Theorem 2.11 for convex, inequality constrained problems can be generalized to general convex problems as follows:

Theorem 2.19: KKT Sufficient Conditions for Convex Programs

Let $f : ℝ^{n_{x}} \to ℝ$ , $g_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{g}$ be convex and differentiable functions. Let also $h_{i} : ℝ^{n_{x}} \to ℝ, i = 1, \dots, n_{h}$ be affine functions. Consider the problem to minimize $f (x)$ subject to $x \in S : = {x \in ℝ^{n_{x}} : g (x) \leq 0, h (x) = 0}$ . If $(x ⋆, ν ⋆, λ ⋆)$ satisfies the KKT conditions of Theorem 2.17 then $x ⋆$ is a global minimizer for $f$ on $S$ .

Theorem 2.20: Farkas’ Theorem

Let $A \in ℝ^{m \times n}$ and $c \in ℝ^{n}$ . Then, exactly one of the following two statements holds:

System 1. $\exists x \in ℝ^{n} such that A x \leq 0 and c ⊤ x > 0$ ,
System 2. $\exists y \in ℝ^{n} such that A^{⊤} y = c and y \geq 0$ .

Proof. See, e.g., [6, Theorem 2.4.5] for a proof. □

Farkas’ Theorem is used extensively in the derivation of optimality conditions of (linear and) nonlinear programming problems. A geometrical interpretation of Farkas’ Theorem is shown in Fig. 1.A.1.. If a1, . . . , am denote the rows of A, then system 2 has a solution if c lies in the convex cone generated by a1, . . . , am; On the other hand, system 1 has a solution if the closed convex cone x : Ax 0 and the open half-space x : cTx ¿ 0 have a nonempty intersection.

Corollary 2.1: Gordan’s Theorem

Let $A \in ℝ^{m \times n}$ . Then, exactly one of the following two statements holds:

System 1. $\exists x \in ℝ^{n} such that A x < 0$ ,
System 2. $\exists y \in ℝ^{m} y \neq 0 such that A^{⊤} y = 0 and y \geq 0$ .

Proof. System 1 can be equivalently written as $A x + ρ e$ where $ρ > 0$ is a scalar and $e$ is a vector of $m$ ones. Rewriting this in the form of System 1 in Farkas’ Theorem 2.20, we get $(A e) p$ and $(0, . . ., 0, 1) p > 0$ where $p : = (x ρ)$ . The associated System 2 by Theorem 2.20 states that ${(A e)}^{⊤} {(0, . . ., 0, 1)}^{⊤}$ and $y \geq 0$ for some $y \in ℝ^{m}$ , i.e., $A^{⊤} y = 0$ , $e ⊤ y = 1$ and $y \geq 0$ , which is equivalent to the System 2 of the corollary. □

Lemma 2.4

Let $P$ and $Q$ be two symmetric matrices, such that $Q ≽ 0$ and $P ≻ 0$ on the null space of $Q$ (i.e., $y ⊤ P y > 0 \forall y \neq 0 such that Q y = 0$ ). Then,

\exists \bar{c} > 0 such that P + c Q ≻ 0 \forall c > \bar{c}

Proof. Assume the contrary. Then,

\forall k > 0, \exists x k, ∥ x k ∥ = 1 such that x k ⊤ P x k + k x k ⊤ Q x k \leq 0

Consider a subsequence ${{x}_{k}} 𝒦$ converging to some $\bar{x}$ with $∥ \bar{x} ∥ = 1$ . Dividing the last equation by $k$ , and taking the limit as $k \in 𝒦 \to \infty$ , we obtain

{\bar{x}}^{⊤} Q \bar{x} \leq 0

On the other hand, $Q$ being positive semidefinite, we must have

{\bar{x}}^{⊤} Q \bar{x} \geq 0

Hence ${\bar{x}}^{⊤} Q \bar{x} = 0$ . That is, using the hypothesis, ${\bar{x}}^{⊤} P \bar{x} > 0$ . This contradicts the fact that

{\bar{x}}^{⊤} P \bar{x} + \underset{k \to 0 k \in 𝒦}{limsup} k x k ⊤ Q x k \leq 0

□

2.9 Numerical methods for nonlinear programming problems

Nowadays, strong and efficient mathematical programming techniques are available for solving a great variety of nonlinear problems, which are based on solid theoretical results and extensive numerical studies. Approximated functions, derivatives and optimal solutions can be employed together with optimization algorithms to reduce the computational time. The aim of this section is not to describe state-of-the-art algorithms in nonlinear programming, but to explain, in a simple way, a number of modern algorithms for solving nonlinear problems. These techniques are typically iterative in the sense that, given an initial point $x 0$ , a sequence of points, ${x k}$ , is obtained by repeated application of an algorithmic rule. The objective is to make this sequence converge to a point $\bar{x}$ in a pre-specified solution set. Typically, the solution set is specified in terms of optimality conditions, such as those presented in 2.5 through 2.8.

We start by recalling a number of concepts in 2.9.1. Then, we discuss the principles of Newton-like algorithms for nonlinear systems in 2.9.2, and use these concepts for the solution of unconstrained optimization problems in 2.9.3. Finally, algorithms for solving general nonlinear problems are presented in 2.9.7, with emphasizes on sequential unconstrained minimization (SUM) and sequential quadratic programming (SQP) techniques.

2.9.1 Preliminaries

Two essential questions must be addressed concerning iterative algorithms. The first question, which is qualitative in nature, is whether a given algorithm in some sense yields, at least in the limit, a solution to the original problem; the second question, the more quantitative one, is concerned with how fast the algorithm converges to a solution. We elaborate on these concepts in this subsection.

The convergence of an algorithm is said to be asymptotic when the solution is not achieved after a finite number of iterations; except for particular problems such as linear and quadratic programming, this is generally the case in nonlinear programming. That is, a very desirable property of an optimization algorithm is global convergence:

Definition 2.23: Global Convergence, Local Convergence

An algorithm is said to be globally convergent if, for any initial point $x 0$ , it generates a sequence of points that converges to a point $\bar{x}$ in the solution set. It is said to be locally convergent if there exists a $ρ > 0$ such that for any initial point $x 0$ such that $‖ \bar{x} - x 0 ‖ < ρ$ it generates a sequence of points converging to $\bar{x}$ in the solution set.

Most modern mathematical programming algorithms are globally convergent. Locally convergent algorithms are not useful in practice because the neighborhood of convergence is not known in advance and can be arbitrarily small.

Next, what distinguishes optimization algorithms with the global convergence property is the order of convergence:

Definition 2.24: Order of Convergence

The order of convergence of a sequence ${x k} \to \bar{x}$ is the largest nonnegative integer p such that

\lim_{k \to \infty} \frac{∥ x k + 1 - \bar{x} ∥}{∥ x k - \bar{x}} = β < \infty

When $p = 1$ and the convergence ratio $β < 1$ , the convergence is said to be linear; if $β = 0$ , the convergence is said to be superlinear. When $p = 2$ , the convergence is said to be quadratic.

Since they involve the limit when $k \to \infty$ , $p$ and $β$ are a measure of the asymptotic rate of convergence, i.e., as the iterates gets close to the solution; yet, a sequence with a good order of convergence may be very ”slow” far from the solution. Clearly, the convergence is faster when $p$ is larger and $β$ is smaller. Near the solution, if the convergence rate is linear, then the error is multiplied by $β$ at each iteration. The error reduction is squared for quadratic convergence, i.e., each iteration roughly doubles the number of significant digits. The methods that will be studied hereafter have convergence rates varying between linear and quadratic.

Example 2.19. Consider the problem to minimize $f (x) = x^{2}$ , subject to $x \geq 1$ .

Let the (point-to-point) algorithmic map $ℳ_{1}$ be defined defined as $ℳ_{1} (x) = \frac{1}{2} (x + 1)$ . It is easily verified that the sequence obtained by applying the map $ℳ_{1}$ , with any starting point, converges to the optimal solution $x^{⋆} = 1$ , i.e., $ℳ_{1}$ is globally convergent. For example, with $x^{0} = 4$ , the algorithm generates the sequence ${4, 2.5, 1.75, 1.375, 1.1875, . . .}$ . We have $(x^{k + 1} - 1) = \frac{1}{2} (x^{k} - 1)$ , so that the limit in Definition 2.24 is $β = \frac{1}{2}$ with $p = 1$ ; moreover, for $p > 1$ , this limit is infinity. Consequently, $x^{k} \to 1$ linearly with convergence ratio $\frac{1}{2}$ .

On the other hand, consider the (point-to-point) algorithmic map $ℳ_{2}$ be defined defined as $ℳ_{2} (x; k) = 1 + \frac{1}{2^{k + 1}} (x - 1)$ . Again, the sequence obtained by applying $ℳ_{2}$ converges to $x^{⋆} = 1$ , from any starting point. However, we now have $\frac{| x^{k + 1} - 1 |}{| x^{k} - 1 |} = \frac{1}{2^{k}}$ , which approaches $0$ as $k \to \infty$ . Hence, $x^{k} \to 1$ superlinearly in this case. With $x^{0} = 4$ , the algorithm generates the sequence ${4, 2.5, 1.375, 1.046875 \dots}$ . The algorithmic maps $ℳ_{1}$ and $ℳ_{2}$ are illustrated on the left and right plots in Fig 1.14., respectively.

It should also be noted that for most algorithms, the user must set initial values for certain parameters, such as the starting point and the initial step size, as well as parameters for terminating the algorithm. Optimization procedures are often quite sensitive to these parameters, and may produce different results, or even stop prematurely, depending on their values. Therefore, it is crucial for the user to understand the principles of the algorithms used, so that he or she can select adequate values for the parameters and diagnose the reasons of a premature termination (failure).

2.9.2 Newton-like Algorithms for nonlinear Systems

The fundamental approach to most iterative schemes was suggested over 300 years ago by Newton. In fact, Newton’s method is the basis for nearly all the algorithms that are described herein.

Suppose one wants to find the value of the variable $x \in ℝ^{n_{x}}$ such that

ϕ (x) = 0

where $ϕ : ℝ^{n_{x}} \to ℝ^{n_{x}}$ is continuously differentiable. Let us assume that one such solution exists, and denote it by $x ⋆$ . Let also $x$ be a guess for the solution. The basic idea of Newton’s method is to approximate the nonlinear function $ϕ$ by the first two terms in its Taylor series expansion about the current point $x$ . This yields a linear approximation for the vector function $ϕ$ at the new point $\bar{x}$ ,

ϕ (\bar{x}) = ϕ (x) + \frac{\partial ϕ}{\partial x} [\bar{x} - x] = ϕ (x) + \nabla ϕ {(x)}^{⊤} [\bar{x} - x]

Using this linear approximation, and provided that the Jacobian matrix $\frac{\partial ϕ}{\partial x}$ is nonsingular, a new estimate for the solution $x ⋆$ can be computed by solving the previous equation for $ϕ (x) = 0$

\bar{x} = x - (\frac{\partial ϕ}{\partial x} | x)^{- 1} ϕ (x)

Letting $d = - (\frac{\partial ϕ}{\partial x} | x)^{- 1} ϕ (x)$ , we get the update $\bar{x} = x + d$ .

Of course, $ϕ$ being a nonlinear function, one cannot expect that $ϕ (\bar{x}) = 0$ , but there is much hope that $\bar{x}$ be a better estimate for the root $x ⋆$ than the original guess $x$ . In other words, we might expect that

‖ \bar{x} - x ⋆ ‖ \leq ‖ x - x ⋆ ‖ and ‖ ϕ (\bar{x}) ‖ \leq ‖ ϕ (x) ‖

If the new point is an improvement, then it makes sense to repeat the process, thereby defining a sequence of points ${x 0, x 1, \dots}$ . An algorithm implementing Newton’s method is as follows:

Algorithm 2.1

Initialization Step:

Let $𝜖 > 0$ be a termination scalar, and choose an initial point $x 0$ .

Let $k = 0$ and go to the main step.

Main Step:

Solve the linear system $(\frac{\partial ϕ}{\partial x} | x k) d k = - ϕ (x k)$ for $d k$ .
Compute the new estimate $x k + 1 = x k + d k$
If $‖ ϕ (x k + 1) ‖ < 𝜖$ , stop; otherwise, replace $k \leftarrow k + 1$ , and go to step 1.

It can be shown that the rate of convergence for Newton’s method is quadratic (see Definition 2.24). Loosely speaking, it implies that each successive estimate of the solution doubles the number significant digits, which is a very desirable property for an algorithm to possess.

Theorem 2.21: Quadratic convergence for Newton’s algorithm

Let $ϕ : ℝ^{n_{x}} \to ℝ^{n_{x}}$ be continuously differentiable, and consider Newton’s algorithm defined by the map $ℳ (x) : = x - \nabla ϕ {(x)}^{- ⊤} ϕ (x)$ . Let $x ⋆$ be such that $ϕ (x ⋆) = 0$ , and suppose that $\nabla ϕ {(x}^{⋆}) ⊤$ is nonsingular. Let the starting point $x 0$ be sufficiently close to $x ⋆$ , so that there exist $c_{1}, c_{2} > 0$ with $c_{1} c_{2} ‖ x 0 - x ⋆ ‖ < 1$ , and

\begin{aligned} ‖ \nabla ϕ {(x)}^{- T} ‖ \leq c_{1} \end{aligned}

^a ^a For a rectangular matrix A ∈ ℝ n × m we define the operator norm ‖ A ‖ = sup { ‖ A x ‖ , x ∈ ℝ m such that ‖ x ‖ = 1 } ‖ ϕ ( x ⋆ ) − ϕ ( x ) − ∇ ⁡ ϕ ( x ) ⊤ ⁡ [ x ⋆ − x ] ‖ ≤ c 2 ‖ x ⋆ − x ‖ 2

for each $x$ satisfying $‖ x ⋆ - x ‖ \leq ‖ x ⋆ - x 0 ‖$ . Then, Newton’s algorithm converges with a quadratic rate of convergence.

Proof. See [6, Theorem 8.6.5] for a proof. □

But can anything go wrong with Newton’s method?

Lack of Global Convergence First and foremost, if the initial guess is not sufficiently close to the solution, i.e., within the region of convergence, Newton’s method may diverge. Said differently, Newton’s method as presented above does not have the global convergence property (see Definition 2.23 and Example 2.20 hereafter). This is because $d k : = \nabla ϕ {(x}^{k}) - ⊤ ϕ (x k)$ may not be a valid descent direction far from the solution, and even if , a unit step size might not give a descent in $ϕ$ . Globalization strategies, which aim at correcting this latter deficiency, will be presented in 2.9.4 in the context of unconstrained optimization.

Singular Jacobian Matrix A second difficulty occurs when the Jacobian matrix $\nabla ϕ (x k)$ becomes singular during the iteration process, since the correction defined by the Newton’s algorithm is not defined in this case. Note that the assumption in Theorem 2.21 guarantees that $\nabla ϕ {(x)}^{⊤}$ cannot be singular. But when the Jacobian matrix is singular at the solution point $x ⋆$ , then Newton’s method looses its quadratic convergence property.

Computational Efficiency Finally, at each iteration, Newton’s method requires (i) that the Jacobian matrix $\nabla ϕ (x k)$ be computed, which may be difficult and/or costly especially when the expression of $ϕ (x)$ is complicated, and (ii) that a linear system be solved. The analytic Jacobian can be replaced by a finite-difference approximation, yet this is costly as $n_{x}$ additional evaluations of $ϕ$ are required at each iteration. With the objective of reducing the computational effort, quasi-Newton methods generate an approximation of the Jacobian matrix, based on the information gathered from the iteration progress. To avoid solving a linear system for the search direction, variants of quasi-Newton methods also exist that generate an approximation of the inverse of the Jacobian matrix. Such methods will be described in 2.9.5 in the context of unconstrained optimization.

Example 2.20. Consider the problem to find a solution to the nonlinear equation

f (x) = \arctan (x) = 0

The Newton iteration sequence obtained by starting from $x^{0} = 1$ is as follows:

\begin{matrix} ̲ & ̲ & ̲ \\ k & x^{k} & | f (x^{k}) | \\ ̲ & ̲ & ̲ \\ 0 & 1 & 0.785398 \\ 1 & - 0.570796 & 0.518669 \\ 2 & 0.116860 & 0.116332 \\ 3 & - 1.061022 \times 1 0^{- 3} & 1.061022 \times 1 0^{- 3} \\ 4 & 7.963096 \times 1 0^{- 10} & 7.963096 \times 1 0^{- 10} \\ ̲ & ̲ & ̲ \end{matrix}

Notice the very fast convergence to the solution $x^{⋆} = 0$ , as could be expected from Theorem 2.21 . The first three iterations are represented in Fig. 1.15., on the left plot.

However, convergence is not guaranteed for any initial guess. More precisely, it can be shown that Newton’s method actually diverges when the initial guess is chosen such that $| x^{0} | > α$ , with $α = 1.3917452002707$ being a solution of $a r c t a n (z) = \frac{2 z}{1 + z^{2}}$ ; further, the method cycles indefinitely for $x^{0} = α$ . Both these situations are illustrated in the right plot and the bottom plot of Fig. 1.15., respectively.

2.9.3 Unconstrained Optimization

We now turn to a description of basic techniques used for iteratively solving unconstrained problems of the form:

min f (x) x \in ℝ^{n_{x}}

Many unconstrained optimization algorithms work along the same lines. Starting from an initial point, a direction of movement is determined according to a fixed rule, and then a move is made in that direction so that the objective function value is reduced; at the new point, a new direction is determined and the process is repeated. The main difference between these algorithms rest with the rule by which successive directions of movement are selected. A distinction is usually made between those algorithms which determine the search direction without using gradient information (gradient-free methods), and those using gradient (and higher-order derivatives) information (gradient-based methods). Here, we shall focus our attention on the latter class of methods, and more specifically on Newton-like algorithms.

Throughout this subsection, we shall assume that the objective function $f$ is twice continuously differentiable. By theorem 2.3, a necessary condition for $x ⋆$ to be a local minimum of $f$ is $\nabla f (x ⋆) = 0$ . Hence, the idea is to devise an iterative scheme that finds a point satisfying the foregoing condition. Following the techniques discussed earlier in 2.9.2, this can be done by using a Newton-like algorithm, with $ϕ$ corresponding to the gradient $\nabla f$ of $f$ , and $\nabla ϕ ⊤$ to its Hessian matrix $\nabla^{2} f$ .

At each iteration, a new iterate $x k + 1$ is obtained such that the linear approximation to the gradient at that point is zero,

\nabla f (x k + 1) \approx \nabla f (x k) + \nabla^{2} f (x k) [x k + 1 - x k] = 0

The linear approximation yields the Newton search direction:

d k : = x k + 1 - x k = - \nabla^{2} f {(x}^{k}) - 1 \nabla f (x k)

As discussed in 2.9.2, if it converges, Newton’s method exhibits a quadratic rate of convergence when the Hessian matrix $\nabla^{2} f (x ⋆)$ is nonsingular at the solution point. However, since the Newton iteration is based on finding a zero of the gradient vector, there is no guarantee that the step will move towards a local minimum, rather than a saddle point or even a maximum. To preclude this, the Newton steps should be taken downhill, i.e., the following descent condition should be satisfied at each iteration,

\nabla f (x k) ⊤ d k < 0

Interestingly enough, with the Newton direction iteration, the descent condition becomes

- \nabla f {(x}^{k}) ⊤ \nabla^{2} f {(x}^{k}) - 1 \nabla f (x k) < 0

That is, a sufficient condition to obtain a descent direction at $x k$ is that the Hessian matrix $\nabla^{2} f (x k)$ be positive definite. Moreover, if $\nabla^{2} f (x ⋆)$ is positive definite at a local minimizer $x ⋆$ of $f$ , then the Newton iteration converges to $x ⋆$ when started sufficiently close to $x ⋆$ . (Recall that, by Theorem 2.8, positive definiteness of $\nabla^{2} f (x ⋆)$ is a sufficient condition for a local minimum of $f$ to be a strict local minimum.)

We now discuss two important improvements to Newton’s method, which are directly related to the issues discussed in Subsection 2.9.2, namely (i) the lack of global convergence, and (ii) computational efficiency.

2.9.4 Globalization Strategies

Up to this point, the development has focused on the application of Newton’s method. However, even in the simplest one-dimensional applications, Newton’s method has deficiencies (see, e.g., Example 2.20). Methods for correcting global convergence deficiencies are referred to as globalization strategies. It should be stressed than an efficient globalization strategy should only alter the iterates when a problem is encountered, but it should not impede the ultimate behavior of the method, i.e., the quadratic convergence of a Newton’s method should be retained.

In unconstrained optimization, one can detect problems in a very simple fashion, by monitoring whether the next iterate $x k + 1$ satisfies a descent condition with respect to the actual iterate $x k$ , e.g., $f (x k + 1) < f (x k)$ . Then, either one of two globalization strategies can be used to correct the Newton step. The first strategy, known as line search method, is to alter the magnitude of the step; the second one, known as trust region method, is to modify both the step magnitude and direction. We shall only concentrate on the former class of globalization strategies subsequently.

A line search method proceeds by replacing the full Newton step $x k + 1 = x k + α d k$ with

x k + 1 = x k + α d k

where the step-length $α \geq 0$ is chosen such that the objective function is reduced,

f (x k + α d k) < f (x k)

Clearly, the optimal choice for $α$ would be $α^{⋆} = argmin {l (α) = f (x k + α d k α \geq 0}$ . The resulting minimization problem can be solved by any one-dimensional exact minimization technique (e.g., Newton’s method itself). However, such techniques are costly in the sense that they often require many iterations to converge and, therefore, many function (or even gradient) evaluations.

In response to this, most modern algorithms implement so-called inexact line search criteria, which aim to find a step-length giving an ”acceptable” decrease in the objective function. Note that sacrificing accuracy, we might impair the convergence of the overall algorithm that iteratively employs such a line search. However, by adopting a line search that guarantees a sufficient degree of descent in the objective function, the convergence of the overall algorithm can still be established.

We now describe one popular definition of an acceptable step-length known as Armijo’s rule. Other popular approaches are the quadratic and cubic fit techniques, as well as Wolfe’s and Glodstein’s tests. Armijo’s rule is driven by two parameters $0 < κ_{1} < 1$ and $κ_{2} > 1$ , which respectively manage the acceptable step-length from being too large or too small (Typical vales are $κ_{1} = 0.2$ and $κ_{2} = 2$ ). Define the line search function $l (α) : = f (x k + α d k)$ , for $α \geq 0$ , and consider the modified first-order approximation $\hat{l} (α) = l (0) + κ_{1} α l^{'} (0)$ . A step-length $\bar{α} \in (0, 1)$ is deemed acceptable if the following conditions hold:

\begin{aligned} l (\bar{α}) \leq \bar{l} (\bar{α}) \\ l (κ_{2} \bar{α}) \geq \hat{l} (κ_{2} \bar{α}) \end{aligned}

The first condition prevents the step-length $\bar{α}$ from being too large, whereas the second condition prevents $\bar{α}$ from being too small. The acceptable region defined by the Armijo’s rule is shown in Fig. 1.16. below.

2.9.5 Recursive Updates

Another limitation of Newton’s method when applied to unconstrained optimization problems is that the Hessian matrix of the objective function is needed at each iteration, then a linear system must be solved for obtaining the search direction. For many applications, this can be a costly computational burden. In response to this, quasi-Newton methods attempt to construct this information recursively. However, by so doing, the quadratic rate of convergence is lost.

The basic idea for many quasi-Newton methods is that two successive iterates $x k, x k + 1$ together with the corresponding gradients $\nabla f (x k), \nabla f (x k + 1)$ , yield curvature information by means of the first-order approximation relation

\nabla f (x k + 1) = \nabla f (x k) + \nabla^{2} f (x k) d k + o (‖ d ‖)

with $d k = x k + 1 - x k$ . In particular, given $n_{x}$ linearly independent iteration increments $δ 0, \dots, δ n_{x} - 1$ an approximation of the Hessian matrix can be obtained as

\nabla^{2} f (x n_{x}) \approx [γ 0 \dots γ n_{x} - 1] {[d}^{0} \dots d n_{x} - 1] - 1

or for the inverse Hessian matrix as

\nabla^{2} f (x n_{x}) \approx [d 0 \dots d n_{x} - 1] {[γ}^{0} \dots γ n_{x} - 1] - 1

where $γ k = \nabla f (x k + 1) - \nabla f (x k)$

Note that when the objective function is quadratic, the previous relations are exact. Many interesting quasi-Newton methods use similar ways, although more sophisticated, to construct an approximate Hessian matrix $B^{k}$ that progressively approaches the inverse Hessian. One of the most popular class of quasi-Newton methods (known as the Broyden family) proceeds as follows:

\begin{aligned} B^{k + 1} & : = B^{k} + \frac{d}{k} - \frac{B^{k} γ}{k} + \\ + ξ γ k^{⊤} B^{k} γ k (\frac{d}{k} d k^{⊤} γ k - \frac{B^{k} γ}{k}) (\frac{d}{k} d k^{⊤} γ k - \frac{B^{k} γ}{k}) ⊤ \end{aligned}

where $0 \leq ξ \leq 1$ . It is easily seen that when supplemented with a line search strategy, $d k^{⊤} γ k < 0$ at each $k$ , and hence the Hessian matrix approximations are guaranteed to exist. Moreover, it can be shown that the successive approximates remain positive definite provided that $B^{0}$ is itself positive definite.

By setting $ξ = 0$ , the previous equation yields the Davidon-Fletcher-Powell (DFP) method, which is historically the first quasi-Newton method, while setting $ξ = 1$ gives the Broyden-Fletcher-Goldfard-Shanno (BFGS) method, for which there is substantial evidence that it is the best general purpose quasi-Newton method currently known.

2.9.6 Summary

A Newton-like algorithm including both a line search method (Armijo’s rule) and Hessian recursive update (DFP update) is as follows:

Algorithm 2.2

Initialization Step:

Let $𝜖 > 0$ be a termination scalar, and choose an initial point $x 0 \in R^{n_{x}}$ and a symmetric, positive definite matrix $B^{0} \in ℝ^{n_{x} \times n_{x}}$ . Let $k = 0$ , and go to the main step.

Main Step:

Search Direction - $d k = - B^{k} \nabla f (x k)$ .
Line Search - Find a step $α^{k}$ satisfying Armijo’s conditions.

Update - Compute the new estimates:

\begin{aligned} x k + 1 : = x k + α^{k} d k \\ B^{k + 1} : = B^{k} + \frac{d}{k} - \frac{B^{k} γ}{k} \end{aligned}

with $d k : = x k + 1 - x k$ and $γ k : = \nabla f (x k + 1) - \nabla f (x k)$

If $‖ \nabla f (x k + 1 ‖ < 𝜖$ , stop; otherwise, replace $k \leftarrow k + 1$ , and go to step 1.

The standard unconstrained optimization algorithm in MATLAB (i.e. the function fminunc) is an implementation of quasi-Newton’s method, with DFP or BFGS update, and a line search strategy.

Example 2.21. Consider the problem to find a minimum to Rosenbrock’s function

f (x) = {(1 - x_{1})}^{2} + c {(x_{2} - x_{1}^{2})}^{2}

for $x \in ℝ^{2}$ with $c : = 105$ . We solved this problem using the function fminunc .

The results are shown in Fig. 1.17. Observe the slow convergence of the iterates far from the optimal solution $x ⋆ = {(1, 1)}^{⊤}$ but the very fast convergence in the vicinity of $x ⋆$ .

2.9.7 Constrained Nonlinear Optimization

In this subsection, we turn our attention to algorithms for iteratively solving constrained problems of the form:

min f (x); x \in X \subset ℝ^{n_{x}}

Many modern deterministic algorithms for constrained NLP problems are based on the (rather natural) principle that, instead of solving a difficult problem directly, one had better solve a sequence of simpler, but related, subproblems, which converges to a solution of the original problem either in a finite number of steps or in the limit. Working along these lines, two classes of algorithms can be distinguished for solution of NLP problems with equality and/or inequality constraints. On the one hand, penalty function and interior-point methods consist of solving the problem as a sequence of unconstrained problems (or problems with simple constraints), so that algorithms for unconstrained optimization can be used. These methods, which do not rely on the KKT theory described earlier in 2.6 through 2.8, shall be briefly presented in 2.9.8 and 2.9.9. On the other hand, Newton-like methods solve NLP problems by attempting to find a point satisfying the necessary conditions of optimality (KKT conditions in general). Sequential quadratic programming (SQP), which shall be presented in 2.9.10, represents one such class of methods.

2.9.8 Penalty Function Methods

Methods using penalty functions transform a constrained problem into a single unconstrained problem or a sequence of unconstrained problems. This is done by placing the constraints into the objective function via a penalty parameter in a way that penalizes any violation of the constraints. To illustrate it, consider the NLP problem:

\begin{matrix} \min & f (x) \\ s.t. & g (x) \leq 0 \\ h (x) = 0 \\ x \in X \end{matrix}

where $X \subset ℝ^{n_{x}}$ , $x \in ℝ^{n_{x}}$ and $f : X \to ℝ$ , $g : X \to ℝ^{n_{g}}$ and $h : X \to ℝ^{n_{h}}$ are defined on $X$ .

In general, a suitable penalty function $α (x)$ for the minimization problem is defined by:

α (x) = \sum_{k = 1}^{n_{g}} ϕ [g_{k} (x)] + \sum_{k = 1}^{n_{h}} ψ [h_{k} (x)]

where $ϕ$ and $ψ$ are continuous functions satisfying the conditions:

{\begin{matrix} ϕ (z) = 0 & if z \leq 0 \\ ϕ (z) > 0 & otherwise \end{matrix} and {\begin{matrix} ψ (z) = 0 & if z = 0 \\ ψ (z) > 0 & otherwise \end{matrix}

Typically, $ϕ$ and $ψ$ are of the forms

ϕ (z) = {(\max {0, z})}^{p} and ψ (z) = | z |^{p}

with $p$ a positive integer (taking $p \geq 2$ provides continuously differentiable penalty functions). The function $f (x) + μ α (x)$ is referred to as the auxiliary function.

Example 2.22. Consider the problem to minimize $f (x) = x$ , subject to $g (x) = - x + 2 \leq 0$ . It is immediately evident that the optimal solution lies at the point $x^{⋆} = 2$ , and has objective value $f (x^{⋆}) = 2$ .

Now, consider the penalty problem to minimize $f (x) + μ α (x) = x + μ \max {0, 2 - x}^{2}$ in $ℝ$ , where $μ$ is a large number. Note first that for any $μ$ , the auxiliary function is convex since it is sum of convex functions. Thus, a necessary and sufficient condition for optimality is that the gradient of $f (x) + μ α (x)$ be equal to zero, yielding $x^{μ} = 2 - \frac{1}{2 μ}$ . Thus, the solution of the penalty problem can be made arbitrarily close to the solution of the original problem by choosing $μ$ sufficiently large. Moreover, $f (x^{μ}) + μ α (x^{μ}) = 2 - \frac{1}{4 μ}$ , which can also be made arbitrarily close to $f (x^{⋆})$ by taking $μ$ sufficiently large. These considerations are illustrated in Fig. 1.18. below.

The conclusions of Example 2.22 that the solution of the penalty problem can be made arbitrarily close to the solution of the original problem, and the optimal auxiliary function value arbitrarily close to the optimal objective value, by choosing sufficiently large, is formalized in the following:

Theorem 2.22

Consider a general NLP problem, where $f$ , $g$ and $h$ are continuous functions on $ℝ^{n_{x}}$ and $X$ is a nonempty set in $ℝ^{n_{x}}$ . Suppose that the NLP problem has a feasible solution, and let $α$ be a continuous penalty function. Suppose further that for each $μ$ there exists a solution $x μ \in X$ to the problem $\min {f (x) + μ α (x) : x \in X}$ , and that ${x μ}$ is contained in a compact subset of $X$ . Then,

\min {f (x) : g (x) \leq 0, h (x) = 0, x \in X} = \sup_{μ \geq 0} 𝜃 (μ) = \lim_{μ \to \infty} 𝜃 (μ)

with $𝜃 (μ) : = f (x μ) + μ α (x μ)$ . Furthermore, the limit $\bar{x}$ of any convergent subsequence of ${x μ}$ is an optimal solution to the original problem and $μ α (x μ) \to 0$ as $μ \to \infty$ .

Proof. See [6, Theorem 9.2.2] for a proof. □

Note that the assumption that $X$ is compact is necessary, for it possible that the optimal objective values of the original and penalty problems are not equal otherwise. Yet, this assumption is not very restrictive in most practical cases as the variables usually lie between finite upper and lower bounds. Note also that no restriction is imposed on $f$ , $g$ and $h$ other than continuity. However, the application of an efficient solution procedure for the (unconstrained) auxiliary problems may impose additional restriction on these functions (see 2.9.3).

Under the conditions that (i) $f$ , $g$ , $h$ and $ϕ$ , $ψ$ are continuously differentiable, and (ii) $\bar{x}$ is a regular point the solution to the penalty problem can be used to recover the Lagrange multipliers associated with the constraints at optimality. In the particular case where $X = ℝ^{n_{x}}$ , we get

\begin{aligned} ν_{i}^{μ} & = μ ϕ^{'} [g_{i} (x μ)] \forall i \in 𝒜 (\bar{x}) \\ λ_{i}^{μ} & = μ ψ^{'} [h_{i} (x μ)] \forall i = 1, \dots, n_{h} \end{aligned}

The larger $μ$ , the better the approximation of the Lagrange multipliers:

ν μ \to ν ⋆ and λ μ \to λ ⋆ as μ \to \infty

Example 2.23. Consider the same problem as in Example 1.71. The auxiliary function $f (x) + μ α (x) = x + μ \max {0, 2 - x}^{2}$ being continuously differentiable, the Lagrange multiplier associated to the inequality constraint $g (x) = - x + 2 \leq 0$ can be recovered as $ν^{μ} = 2 μ \max {0, 2 - x}$ (assuming $μ > 0$ ). Note that the exact value of the Lagrange multiplier is obtained for each $μ > 0$ here, because $g$ is a linear constraint.

From a computational viewpoint, superlinear convergence rates might be achievable, in principle, by applying Newton’s method (or its variants such as quasi-Newton methods). Yet, one can expect ill-conditioning problems when $μ$ is taken very large in the penalty problem. With a large $μ$ , more emphasis is placed on feasibility, and most procedures for unconstrained optimization will move quickly towards a feasible point. Even though this point may be far from the optimum, both slow convergence and premature termination can occur due to very small step size and finite precision computations (round-off errors).

As a result of the above mentioned difficulties associated with large penalty parameters, most algorithms using penalty functions employ a sequence of increasing penalty parameters. With each new value of the penalty parameter, an optimization technique is employed, starting with the optimal solution corresponding to the previously chosen parameters value. Such an approach is often referred to as sequential unconstrained minimization (SUM) technique. This way, a sequence of infeasible points is typically generated, whose limit is an optimal solution to the original problem (hence the term exterior penalty function approach).

To conclude our discussion on the penalty function approach, we give an algorithm to solve a general NLP problem, where the penalty function used are of the form specified by $ϕ$ and $ψ$

Algorithm 2.3

Initialization Step:

Let $𝜖 > 0$ be a termination scalar, and choose an initial point $x 0$ , a penalty parameter $μ^{0} > 0$ , and a scalar $β > 1$ . Let $k = 0$ and go to the main step.

Main Step:

Starting with

x k

, get a solution to the problem:

x k + 1 \in \arg \min {f (x) + μ^{k} α (x) : x \in X}

If $μ^{k} α (x k + 1) < 𝜖$ , stop; otherwise, let $μ^{k + 1} = β μ^{k}$ , replace $k \leftarrow k + 1$ , and go to step 1.

2.9.9 Interior-Point Methods

Similar to penalty functions, barrier functions can also be used to transform a constrained problem into an unconstrained problem (or into a sequence of unconstrained problems). These functions act as a barrier and prevent the iterates from leaving the feasible region. If the optimal solution occurs at the boundary of the feasible domain, the procedure moves from the interior to the boundary of the domain, hence the name interior-point methods. To illustrate these methods, consider the NLP problem:

\begin{matrix} \min_{x} & f (x) \\ s.t. & g (x) \leq 0 \\ x \in X \end{matrix}

where $X$ is a subset of $ℝ^{n_{x}}$ , and $f : X \to ℝ$ , $g : X \to ℝ^{n_{g}}$ are continuous on $ℝ^{n_{x}}$ . Note that equality constraints, if any, should be accommodated within the set $X$ . (In the case of affine equality constraints, one can possibly eliminate them after solving for some variables in terms of the others, thereby reducing the dimension of the problem.) The reason why this treatment is necessary is because barrier function methods require the set ${x \in ℝ^{n_{x}} : g (x) < 0}$ to be nonempty; this would obviously be not possible if the equality constraints $h (x) = 0$ were accommodated within the set of inequalities as $h (x) \leq 0$ and $h (x) \geq 0$ .

A barrier problem formulates as:

\begin{matrix} \inf_{μ} & 𝜃 (μ) \\ s.t. & μ > 0 \end{matrix}

where $𝜃 (μ) : = \inf {f (x) + μ b (x) : g (x) < 0, x \in X}$ . Ideally, the barrier function $b$ should take value zero on the region ${x : g (x) \leq 0}$ , and value $\infty$ on its boundary. This would guarantee that the iterates do not leave the domain ${x : g (x) \leq 0}$ provided the minimization problem started at an interior point. However, this discontinuity poses serious difficulties for any computational procedure. Therefore, this ideal construction of $b$ is replaced by the more realistic requirement that $b$ be nonnegative and continuous over the region ${x : g (x) \leq 0}$ and approach infinity as the boundary is approached from the interior:

b (x) = \sum_{k = 1}^{n_{g}} ϕ [g_{k} (x)]

where $ϕ$ is a continuous function over ${z : z < 0}$ that satisfies the conditions

{\begin{matrix} ϕ (z) \geq 0 if z < 0 \\ \lim_{z \to 0^{-}} ϕ (z) = + \infty \end{matrix}

In particular, $μ b$ approaches the ideal barrier function described above as $μ$ approaches zero.

Typically barrier functions are

b (x) = - \sum_{k = 1}^{n_{g}} \frac{1}{g_{k} (x)} or b (x) = - \sum_{k = 1}^{n_{g}} \ln [\min {1, - g_{k} (x)}]

The following barrier function, known as Frisch’s logarithmic barrier function, is also widely used

b (x) = - \sum_{k = 1}^{n_{g}} \ln [- g_{k} (x)]

Actually the Frisch’s logarithmic barrier $ϕ (z) = - \ln (- z)$ does not satisfy the nonnegativity requirement for $z < - 1$ . However the requirement on $ϕ$ can be relaxed and it is sufficient that $ϕ$ be positive close to $z = 0$ .

The function $f (x) + μ b (x)$ is referred to as the auxiliary function.

Given $μ > 0$ , evaluating $𝜃 (μ) = inf {f (x) + μ b (x) : g (x) < 0, x \in X}$ seems no simpler than solving the original problem because of the constraint $g (x) < 0$ . However, starting the optimization from a point in the region $S : = {x : g (x) < 0} \cap X$ yields an optimal point in $S$ , even when the constraint $g (x) < 0$ is ignored. This is because $b$ approaches infinity as the iterates approach the boundary of ${x : g (x) \leq 0}$ from within $S$ , hence preventing them from leaving the set $S$ . This is formalized in the following:

Theorem 2.23

Consider a NLP problem with inequality constraints , where $f$ and $g$ are continuous functions on $ℝ^{n_{x}}$ and $X$ is a nonempty closed set in $ℝ^{n_{x}}$ . Suppose that the minimization problem has an optimal solution $x ⋆$ with the property that, given any neighborhood $ℬ_{η} (x ⋆)$ around $x ⋆$ , there exists an $x \in X \cap ℬ_{η} (x ⋆)$ such that $g (x) < 0$ . Suppose further that for each $μ$ , there exists a solution $x μ \in X$ to the problem $\min {f (x) + μ b (x) : x \in X}$ . Then,

\min {f (x) : g (x) \leq 0, x \in X} = \lim_{μ \to 0^{+}} 𝜃 (μ) = \inf_{μ > 0} 𝜃 (μ)

with $𝜃 (μ) : = f (x μ) + μ b (x μ)$ . Furthermore, the limit of any convergent subsequence of ${x μ}$ is an optimal solution to the original problem, and $μ b (x μ) \to 0$ as $μ \to 0^{+}$ .

Proof. See [6, Theorem 9.4.3] for a proof. □

Under the conditions that (i) $f$ , $g$ and $ϕ$ are continuously differentiable, and (ii) $x ⋆$ is a regular point, the solution to the barrier function problem can be used to recover the Lagrange multipliers associated with the constraints at optimality. In the particular case where $X = ℝ^{n_{x}}$ , we get:

ν_{i}^{μ} = μ ϕ^{'} [g_{i} (x μ)] \forall i \in 𝒜 (x ⋆)

The approximation of the Lagrange multipliers, gets better as $μ$ gets closer to $0$ ,

ν μ \to ν ⋆ as μ \to 0^{+}

Example 2.24. Consider the problem to minimize $f (x) = x$ , subject to $g (x) = - x + 2 \leq 0$ , the solution of which lies at the point $x^{⋆} = 2$ with objective value $f (x^{⋆}) = 2$ . Now, consider the barrier function problem to minimize $f (x) + μ b (x) = x - \frac{μ}{2 - x}$ in $ℝ$ , where $μ$ is a large number. Note first that for any $μ$ , the auxiliary function is convex ¹ ¹ $\frac{μ}{2 - x}$ is a concave function on the convex set $S : = {x \in ℝ, x \geq 2}$ , therefore $- \frac{μ}{2 - x}$ is convex and thus $f (x) + μ b (x)$ is sum of convex functions. . Thus, a necessary and sufficient condition for optimality is that the gradient of $f (x) + μ b (x)$ be equal to zero, yielding $x^{μ} = 2 + \sqrt{μ}$ (assuming $μ > 0$ ). Thus, the solution of the penalty problem can be made arbitrarily close to the solution of the original problem by choosing $μ$ sufficiently close to zero. Moreover, $f (x^{μ}) + μ b (x^{μ}) = 2 + 2 \sqrt{μ}$ , which can also be made arbitrarily close to $f (x^{⋆})$ by taking $μ$ sufficiently close to zero. These considerations are illustrated in Fig. 1.19. below. Regarding the Lagrange multiplier associated to the inequality constraint $g (x) = - x + 2 \leq 0$ , the objective and constraint functions being continuously differentiable, an approximate value can be obtained as $ν^{μ} = \frac{μ}{{(2 - x^{μ})}^{2}} = 1$ . Here again, the exact value of the Lagrange multiplier is obtained for each $μ > 0$ because $g$ is a linear constraint.

The use of barrier functions for solving constrained NLP problems also faces several computational difficulties. First, the search must start with a point $x \in X$ such that $g (x) < 0$ , and finding such a point may not be an easy task for some problems. Also, because of the structure of the barrier function $b$ , and for small values of the parameter $μ$ , most search techniques may face serious ill-conditioning and difficulties with round-off errors while solving the problem to minimize $f (x) + μ b (x)$ over $x \in X$ , especially as the boundary of the region ${x : g (x) \leq 0}$ is approached. Accordingly, interior-point algorithms employ a sequence of decreasing penalty parameters ${μ^{k}} \to 0$ as $k \to \infty$ ; with each new value $μ^{k}$ , an optimal solution to the barrier problem is sought by starting from the previous optimal solution. As in the exterior penalty function approach, it is highly recommended to use suitable second-order Newton or quasi-Newton methods for solving the successive barrier problems.

We describe below a scheme using barrier functions for optimizing a NLP problem.

Algorithm 2.4

Initialization Step:

Let $𝜖 > 0$ be a termination scalar, and choose an initial point $x 0$ , with $g (x 0) < 0$ . Let $μ^{0} > 0$ , $β \in (0, 1)$ , $k = 0$ , and go to the main step.

Main Step:

Starting with

x k

, get a solution to the problem:

x k + 1 \in \arg \min {f (x) + μ^{k} b (x) : x \in X}

If $μ^{k} b (x k + 1) < 𝜖$ , stop; otherwise, let $μ^{k + 1} = β μ^{k}$ , replace $k \leftarrow k + 1$ , and go to step 1.

Note that although the constraint $g (x) < 0$ may be ignored, it is considered in the problem formulation as most line search methods use discrete steps, and a step could lead to a point outside the feasible region (where the value of the barrier function is a large negative number), when close to the boundary. Therefore, the problem can effectively be treated as an unconstrained optimization problem only if an explicit check for feasibility is made.

In recent years, there has been much excitement because some variants of the interior-point algorithm can be shown to be polynomial in time for many classes of convex programs. Moreover, interior-point codes are now proving to be highly competitive with codes based on other algorithms, such as SQP algorithms presented subsequently. A number of free and commercial interior-point solvers are given in Tab. 1.1. below.

2.9.10 Sequential Quadratic Programming

Sequential quadratic programming (SQP) methods, also known as successive, or recursive, quadratic programming, employ Newton’s method (or quasi-Newton methods) to directly solve the KKT conditions for the original problem. As a result, the accompanying subproblem turns out to be the minimization of a quadratic approximation to the Lagrangian function subject to a linear approximation to the constraints. Hence, this type of process is also known as a projected Lagrangian, or the Newton-Lagrange, approach. By its nature, this method produces both primal and dual (Lagrange multiplier) solutions.

To present the concept of SQP, consider first the equality constrained nonlinear program:

\begin{matrix} \min_{x} & f (x) \\ s.t. & h (x) = 0 \end{matrix}

where $x \in ℝ^{n_{x}}$ , and $f : ℝ^{n_{x}} \to ℝ$ , $h : ℝ^{n_{x}} \to ℝ^{n_{h}}$ are twice continuously differentiable. We shall also assume throughout that the equality constraints are linearly independent at a solution $x ⋆$ . (The extension for including inequality constraints is considered subsequently.)

By Theorem 2.13, the first-order necessary conditions of optimality for an equality constrained NLP require a primal solution $x ⋆ \in ℝ^{n_{x}}$ and a Lagrange multiplier vector $λ ⋆ \in ℝ^{n_{h}}$ such that:

\begin{aligned} 0 = \nabla_{x} ℒ (x ⋆, λ ⋆) = \nabla f (x ⋆) + \nabla h (x ⋆) λ ⋆ \\ 0 = \nabla_{λ} ℒ (x ⋆, λ ⋆) = h (x ⋆) = 0 \end{aligned}

where $ℒ (x, λ) : = f (x) + λ ⊤ h (x)$ . Now, consider a Newton-like method to solve the latter system of $n_{x} + n_{h}$ nonlinear equations. Given an iterate $(x k, λ k)$ , a new iterate $(x k + 1, λ k + 1)$ is obtained by solving the first-order approximation:

0 = (\begin{matrix} \nabla \end{matrix} x ℒ (x k, λ k) \nabla λ ℒ (x k, λ k)) + (\begin{matrix} \nabla_{x x}^{2} ℒ (x k, λ k \end{matrix}) \nabla h (x k) \nabla h (x k) ⊤ & 0) (\begin{matrix} x \end{matrix} k + 1 - x k λ k + 1 - λ k)

Substituting $\nabla_{x} ℒ (x k, λ k) = \nabla f (x k) + \nabla h (x k) λ k$ , the first row of the system of equations is:

\nabla_{x x}^{2} ℒ (x k, λ k) [x k + 1 - x k] + \nabla h (x k) [λ k + 1 - λ k] + \nabla f (x k) + \nabla h (x k) λ k = 0

The Lagrange multiplier $λ k$ can be simplified and thus:

\nabla_{x x}^{2} ℒ (x k, λ k) [x k + 1 - x k] + \nabla h (x k) λ k + 1 + \nabla f (x k) = 0

The second row of the system of equations is independent of the Lagrange multipliers, it is:

\nabla h {(x}^{k}) ⊤ [x k + 1 - x k] + h (x k) = 0

Therefore the system can be simplified as:

(\begin{matrix} \nabla_{x x}^{2} ℒ (x k, λ k \end{matrix}) \nabla h (x k) \nabla h (x k) ⊤ & 0) (\begin{matrix} d \end{matrix} k λ k + 1) = - (\begin{matrix} \nabla f (x k) \end{matrix} h (x k))

Where we have defined the Newton step on $x$ as $d k = x k + 1 - x k$ . The simplified system can be solved for $(d k, λ k)$ if a solution exists. Setting $x k + 1 : = x k + d k$ , and incrementing $k$ by $1$ , we can then repeat the process until $d k = 0$ happens to solve the linear system. When this occurs, if at all, we shall have found a stationary point to the equality constrained NLP.

Interestingly enough, a quadratic programming (QP) minimization subproblem can be employed in lieu of the foregoing linear system to find any optimal solution for the linearized NLP problem. Given $(x k, λ k)$ we define the subproblem $QP (x k, λ k)$ as:

\begin{matrix} \min_{d} f (x k) + \nabla f (x k) ⊤ d k + \frac{1}{2} d k^{⊤} \nabla_{x x}^{2} ℒ (x k, λ k) d k \\ s.t. h (x k) + \nabla h (x k) ⊤ d k = 0 \end{matrix}

First, note that an optimum $d k$ to $QP (x k, λ k)$ , if it exists, together with the set of Lagrange multipliers $λ k + 1$ associated with the linearized constraints satisfies the linearized KKT conditions. Second, the objective function of $QP (x k, λ k)$ represents not just a quadratic approximation for $f (x)$ but also incorporates an additional term $\frac{1}{2} d k^{⊤} {\sum_{i = 1}^{n_{h}} λ_{i}^{k} \nabla^{2} h_{i} (x k)} d k$ to represent the curvature of the constraints.

Algorithm 2.5: SQP

Initialization Step:

Choose and initial primal/dual point $(x 0, λ 0)$ , let $k = 0$ , and go to the main step.

Main Step:

Solve the quadratic subproblem $QP (x k, λ k)$ to obtain a solution $d k$ long with a set of Lagrange multipliers $λ k + 1$ .
If $d k = 0$ , then $(d k, λ k + 1)$ satisfies the stationarity conditions of the original NLP problem; stop. Otherwise, let $x k + 1 : = x k + d k$ , replace $k \leftarrow k + 1$ , and go to step 1.

In the case $x ⋆$ is a regular stationary solution for the NLP problem which, together with a set of Lagrange multipliers $λ ⋆$ , satisfies the second-order sufficiency conditions of Theorem 2.16, then the matrix

W : = (\begin{matrix} \nabla_{x x}^{2} ℒ (x ⋆, λ ⋆ \end{matrix}) \nabla h (x ⋆) \nabla h (x ⋆) ⊤ & 0)

can be shown to be nonsingular. Hence, the above rudimentary SQP algorithm exhibits a quadratic rate of convergence by Theorem 2.21.

We now consider the inclusion of inequality constraints $g_{i} (x) \leq 0, i = 1, \dots, n_{h}$ in the NLP problem, that is:

\begin{matrix} \min_{x} & f (x) \\ s.t. & h (x) = 0 \\ g (x) \leq 0 \end{matrix}

where the functions $g_{i}$ are twice continuously differentiable.

Given an iterate $(x k, λ k, ν k)$ , where $λ k$ and $ν k \geq 0$ are the Lagrange multiplier estimates for the equality and inequality constraints, respectively, consider the following QP subproblem as a direct extension of $QP (x k, λ k)$ :

\begin{matrix} \min_{d} & f (x k) + \nabla f (x k) ⊤ d k + \frac{1}{2} d k^{⊤} \nabla_{x x}^{2} ℒ (x k, λ k, ν k) d k \\ s.t. & g (x k) + \nabla g (x k) ⊤ d k \leq 0 \\ h (x k) + \nabla h (x k) ⊤ d k = 0 \end{matrix}

where $ℒ (x, λ, ν) : = f (x) + ν ⊤ g (x) + λ ⊤ h (x)$ . Note that the KKT conditions for $QP (x k, λ k, ν k)$ require that, in addition to primal feasibility, Lagrange multipliers $λ k + 1$ , $ν k + 1$ be found such that:

\begin{aligned} \nabla f (x k) + \nabla_{x x}^{2} ℒ (x k, λ k, ν k) d k + \nabla g (x k) ν k + 1 + \nabla h (x k) λ k + 1 & = 0 \\ [g (x k) + \nabla g (x k) ⊤ d k] ⊤ ν k + 1 & = 0 \end{aligned}

with $ν k + 1 \geq 0$ and $λ k + 1$ unrestricted in sign. Clearly, if $d k = 0$ , then $d k$ together with $λ k + 1$ , $ν k + 1$ yields a KKT solution to the NLP original problem. Otherwise, we set $x k + 1 : = x k + d k$ as before, increment $k$ by $1$ , and repeat the process. Regarding convergence rate, it can be shown that if $x ⋆$ is a regular KKT solution which, together with $λ ⋆$ , $ν ⋆$ satisfies the second-order sufficient conditions of Theorem 2.18, and if $(x 0, λ 0, ν 0)$ is initialized sufficiently close to $(x ⋆, λ ⋆, ν ⋆)$ , then the foregoing iterative procedure shall exhibit a quadratic convergence rate.

The SQP method, as presented thus far, obviously shares the disadvantages of Newton’s method: (i) it requires second-order derivatives $\nabla_{x x}^{2} ℒ$ to be calculated, which in addition might not be positive definite, and (ii) it lacks the global convergence property.

Regarding second-order derivatives, a quasi-Newton positive definite approximation can be used for

\nabla_{x x}^{2} ℒ

. For example, given a positive definite approximation

B^{k}

for

\nabla_{x x}^{2} ℒ

, the quadratic problem

QP (x k, λ k, ν k)

can be solved with

\nabla_{x x}^{2} ℒ

replaced by

B^{k}

. For example, an approximation of the inverse Hessian matrix can be obtained via a Broyden-like procedure, with

γ k

given by

γ^{k} : = \nabla ℒ (x k + 1, λ k + 1, ν k + 1) - \nabla ℒ (x k, λ k + 1, ν k + 1)

It can be shown that this modification to the rudimentary SQP algorithm, similar to the quasi-Newton modification of Newton’s algorithm, looses the quadratic convergence rate property. Instead, it can be shown that the convergence is superlinear when initialized sufficiently close to a solution $(x ⋆, λ ⋆, ν ⋆)$ that satisfies both regularity and second-order sufficiency conditions. However, this superlinear convergence rate is strongly based on the use of unit step sizes (see point (ii) below).

In order to remedy the global convergence deficiency, a globalization strategy can be used, e.g., a line search procedure. Unlike unconstrained optimization problems, however, the choice of a suitable line search or merit function providing a measure of progress is not obvious in the presence of constraints. Two such popular choices of a line search function are:

The

ℓ_{1}

Merit Function:

ℓ_{1} (x; μ) : = f (x) + μ [\sum_{i = 1}^{n_{h}} | h_{i} (x) | + \sum_{i = 1}^{n_{g}} \max {0, g_{i} (x)}]

which satisfies the important property that $x ⋆$ is a local minimizer of $ℓ_{1}$ , provided $(x ⋆, λ ⋆, ν ⋆)$ satisfies the second-order sufficient conditions (see Theorem 2.18) and the penalty parameter $μ$ is so chosen that $μ > | λ_{i}^{⋆} |$ , $i = 1, \dots, n_{h}$ and $μ > ν_{i}^{⋆}$ , $i = 1, \dots, n_{g}$ . Yet, the $ℓ_{1}$ merit function is not differentiable at those $x$ with either $g_{i} (x) = 0$ or $h_{i} (x) = 0$ , and it can be unbounded below even though $x ⋆$ is a local minimizer.

The

ℓ_{2}

Augmented Lagrangian (ALAG) Merit Function:

ℓ_{2} (x, λ, ν; μ) : = f (x) + \sum_{i = 1}^{n_{h}} λ_{i} h_{i} (x) + \frac{μ}{2} \sum_{i = 1}^{n_{h}} {[h_{i} (x)]}^{2} + \frac{1}{2} \sum_{i = 1}^{n_{g}} ψ_{i} (x, ν; μ)

with $ψ_{i} (x, ν; μ) : = \frac{1}{μ} (\max {0, ν_{i} + μ g_{i} (x)}^{2} - ν_{i}^{2})$ , has similar properties to the $l_{1}$ merit function, provided $μ$ is chosen large enough, and is continuously differentiable (although its Hessian matrix is discontinuous). Yet, for $x ⋆$ to be a (local) minimizer of $ℓ_{2} (x, λ, ν; μ)$ , it is necessary that $λ = λ ⋆$ and $ν = ν ⋆$ .

An SQP algorithm including the modifications discussed in (i) and (ii) is as follows:

Algorithm 2.6: Improved SQP

Initialization Step:

Choose and initial primal/dual point $(x 0, λ 0, ν 0)$ , with $ν 0 \geq 0$ , and a positive definite matrix $B^{0}$ . Let $k = 0$ , and go to the main step.

Main Step:

Solve the quadratic subproblem $QP (x k, λ k, ν k)$ , with
$\nabla_{x x}^{2} ℒ (x k, λ k, ν k)$ replaced by $B^{k}$ , to obtain a direction $d k$ along with a set of Lagrange multipliers $(λ k + 1, ν k + 1)$ .
If $d k = 0$ , then $(d k, λ k + 1, ν k + 1)$ satisfies the KKT conditions of the original NLP problem; stop.
Find $x k + 1 : = x k + α^{k} d k$ , where $α^{k}$ improves $ℓ_{1} (x k + α d k)$ over ${α \in ℝ : α > 0}$ [or any other suitable merit function]. Update $B^{k}$ to a positive definite matrix $B^{k + 1}$ [e.g., according to the quasi-Newton update]. Replace $k \leftarrow k + 1$ , and go to step 1.

What Can Go Wrong?

The material presented up to this point was intended to give the reader an understanding of how SQP methods should work. Things do not go so smoothly in practice though. We now discuss a number of common difficulties that can be encountered, and suggest remedial actions to correct the deficiencies. Because real applications may involve more than a single difficulty, the user must be prepared to correct all problems before obtaining satisfactory performance from an optimization software.

Infeasible Constraints
One of the most common difficulties occurs when the NLP problem has infeasible constraints, i.e., the constraints taken all together have no solution. Applying general-purpose SQP software to such problems typically produces one or more of the following symptoms:
- one of the QP subproblem happen to be infeasible, which occurs when the linearized constraints have no solution;
- many NLP iterations produce very little progress;
- the penalty parameters $μ$ of the merit functions growsvery large;
- the Lagrange multipliers become very large.
Although robust SQP software attempt to diagnose this situation, ultimately the only remedy is to reformulate the NLP.
Rank-Deficient Constraints
In contrast to the previous situation, it is possible that the constraints be consistent, but the Jacobian matrix of the active constraints, at the solution point, be either ill-conditioned or rank deficient. This situation was illustrated in Examples 1.43 and 1.55. The application of general-purpose SQP software is likely to produce the following symptoms:
- many NLP iterations produce very little progress;
- the penalty parameters $μ$ of the merit functions grows very large;
- the Lagrange multipliers become very large;
- the rank deficiency in the Jacobian of the active constraints is detected.
Note that many symptoms of rank-deficient constraints are the same as those of inconsistent constraints. It is therefore quite common to confuse this deficiency with inconsistent constraints. Again, the remedy is to reformulate the problem.
Redundant Constraints
A third type of difficulty occurs when the NLP problem contains redundant constraints. Two types of redundancy may be distinguished. In the first type, some of the constraints are unnecessary to the problem formulation, which typically results in the following symptoms:
- the Lagrange multipliers are close to zero;
- the solver has difficulty in detecting the active set.
In the second type, the redundant constraints give rise to rank deficiency, and the problem then exhibits symptoms similar to the rank-deficient case discussed previously. Obviously, the remedy is to reformulate the problem by eliminating the redundant constraints.
Discontinuities
Perhaps the biggest obstacle encountered in the practical application of SQP methods (as well as many other NLP methods including SUM techniques) is the presence of discontinuous behavior. All of the numerical methods described herein assume continuous and differentiable objective function and constraints. Yet, there are many common examples of discontinuous functions in practice, including: IF tests in codes; absolute value, $\min$ , and $\max$ functions; linear interpolation of data; internal iterations such as root finding; etc.

Regarding SQP methods, the standard QP subproblems are no longer appropriate when discontinuities are present. In fact, the KKT necessary conditions simply do not apply! The most common symptoms of discontinuous functions are:
- the iterates converge slowly or, even, diverge;
- the line search takes very small steps $(α \approx 0)$ ;
- the Hessian matrix becomes badly ill-conditioned.
The remedy consists in reformulating discontinuous problems into smooth problems: for absolute value, min, and max functions, tricks can be used that introduce slack variables and additional constraints; linear data interpolation can be replaced by higher order interpolation schemes that are continuous through second derivatives; internal iterations can also be handled via additional NLP constraints; etc.

Inaccurate Gradient Estimation

Any SQP code requires that the user supply the objective function and constraint values, as well as their gradient (and possibly their Hessian too). In general, the user is proposed the option to calculate the gradients via finite differences, e.g.,

\nabla f (x) \approx \frac{f (x + δ x) - f (x)}{δ x}

However, this may cause the problem to stop prematurely. First of all, the choice of the perturbation vector $δ x$ is highly non trivial. If too large a value clearly provides inaccurate estimates, too small a value may also result in very bad estimates due to finite arithmetic precision computations. Therefore, one must try to find a trade-off between these two extreme situations. The difficulty stems from the fact that a trade-off may not necessarily exist if the requested accuracy for the gradient is too high. In other word, the error made in the finite-difference approximation of a gradient cannot be made as small as desired. Further, the maximum accuracy that can be achieved with finite difference is both problem dependent (e.g., badly-scaled functions are more problematic than well-scaled functions) and machine dependent (e.g., double precision computations provides more accurate estimates than single precision computations). Typical symptoms of inaccurate gradient estimates in an SQP code are:

the iterates converge slowly, and the solver may stop prematurely at a suboptimal point (jamming);
the line search takes very small steps $(α \approx 0)$ .

The situation can be understood as follows. Assume that the gradient estimate is contaminated with noise. Then, instead of computing the true value $\nabla ℒ (x)$ , we get $ℒ (x) + 𝜖$ . But since the iteration seeks a point such that $\nabla ℒ (x) = 0$ , we can expect either a degraded rate of convergence or, worse, no convergence at all, because ultimately the gradient will be dominated by noise.

To avoid these problems, the user should always consider providing the gradients explicitely to the SQP solver, instead of relying on finite-difference estimates. For large-scale problems, this is obviously a time-consuming and error-prone task. In response to this, efficient algorithmic differentiation tools (also called automatic differentiation) have been developed within the last fifteen years. The idea behind it is that, given a piece of program calculating a number of function values (e.g., in fortran77 or C language), an auxiliary program is generated that calculates the derivatives of these functions.

Scaling

Scaling affects everything! Poor scaling can make a good algorithm behave badly. Scaling changes the convergence rate, termination tests, and numerical conditioning.

The most common way of scaling a problem is by introducing scaled variables of the form

{\tilde{x}}_{k} : = u_{k} x_{k} + r_{k}

for $k = 1, . . ., n_{x}$ with $u_{k}$ and $r_{k}$ being scale weights and shifts, respectively. Likewise, the objective function and constraints are commonly scaled using

\begin{aligned} \tilde{f} & : = ω_{0} f \\ {\tilde{g}}_{k} & : = ω_{k} g_{k} \end{aligned}

for $k = 1, . . ., n_{g}$ . The idea is to let the optimization algorithm work with the well-scaled quantities in order to improve performance. However, what well-scaled quantities mean is hard to define, although conventional wisdom suggests the following hints:

normalize the independent variables to have the same range, e.g., $0 \leq {\tilde{x}}_{k} \leq 1$ ;
normalize the dependent functions to have the same magnitude, e.g., $\tilde{f} \approx {\tilde{g}}_{1} \approx \dots \approx {\tilde{g}}_{n_{g}} \approx 1$ ;
normalize the rows and columns of the Jacobian to be of the same magnitude;
scale the dependent functions so that the Lagrange multipliers are close to one, e.g., $| λ_{1} | \approx \dots \approx | λ_{n_{g}} | \approx 1$ ; etc.